LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ slatme()

subroutine slatme ( integer  N,
character  DIST,
integer, dimension( 4 )  ISEED,
real, dimension( * )  D,
integer  MODE,
real  COND,
real  DMAX,
character, dimension( * )  EI,
character  RSIGN,
character  UPPER,
character  SIM,
real, dimension( * )  DS,
integer  MODES,
real  CONDS,
integer  KL,
integer  KU,
real  ANORM,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( * )  WORK,
integer  INFO 
)

SLATME

Purpose:
    SLATME generates random non-symmetric square matrices with
    specified eigenvalues for testing LAPACK programs.

    SLATME operates by applying the following sequence of
    operations:

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

    2. If complex conjugate pairs are desired (MODE=0 and EI(1)='R',
         or MODE=5), certain pairs of adjacent elements of D are
         interpreted as the real and complex parts of a complex
         conjugate pair; A thus becomes block diagonal, with 1x1
         and 2x2 blocks.

    3. If UPPER='T', the upper triangle of A is set to random values
         out of distribution DIST.

    4. If SIM='T', A is multiplied on the left by a random matrix
         X, whose singular values are specified by DS, MODES, and
         CONDS, and on the right by X inverse.

    5. If KL < N-1, the lower bandwidth is reduced to KL using
         Householder transformations.  If KU < N-1, the upper
         bandwidth is reduced to KU.

    6. If ANORM is not negative, the matrix is scaled to have
         maximum-element-norm ANORM.

    (Note: since the matrix cannot be reduced beyond Hessenberg form,
     no packing options are available.)
Parameters
[in]N
          N is INTEGER
           The number of columns (or rows) 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, and for the
           upper triangle (see UPPER).
           '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 SLATME
           to continue the same random number sequence.
           Changed on exit.
[in,out]D
          D is REAL array, dimension ( N )
           This array is used to specify the eigenvalues of A.  If
           MODE=0, then D is assumed to contain the eigenvalues (but
           see the description of EI), otherwise they will be
           computed according to MODE, COND, DMAX, and RSIGN and
           placed in D.
           Modified if MODE is nonzero.
[in]MODE
          MODE is INTEGER
           On entry this describes how the eigenvalues are to
           be specified:
           MODE = 0 means use D (with EI) 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.  Each odd-even pair
                    of elements will be either used as two real
                    eigenvalues or as the real and imaginary part
                    of a complex conjugate pair of eigenvalues;
                    the choice of which is done is random, with
                    50-50 probability, for each pair.
           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 between 1 and 4, D has entries ranging
              from 1 to 1/COND, if between -1 and -4, D has entries
              ranging from 1/COND to 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))).  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]EI
          EI is CHARACTER*1 array, dimension ( N )
           If MODE is 0, and EI(1) is not ' ' (space character),
           this array specifies which elements of D (on input) are
           real eigenvalues and which are the real and imaginary parts
           of a complex conjugate pair of eigenvalues.  The elements
           of EI may then only have the values 'R' and 'I'.  If
           EI(j)='R' and EI(j+1)='I', then the j-th eigenvalue is
           CMPLX( D(j) , D(j+1) ), and the (j+1)-th is the complex
           conjugate thereof.  If EI(j)=EI(j+1)='R', then the j-th
           eigenvalue is D(j) (i.e., real).  EI(1) may not be 'I',
           nor may two adjacent elements of EI both have the value 'I'.
           If MODE is not 0, then EI is ignored.  If MODE is 0 and
           EI(1)=' ', then the eigenvalues will all be real.
           Not modified.
[in]RSIGN
          RSIGN is CHARACTER*1
           If MODE is not 0, 6, or -6, and RSIGN='T', then the
           elements of D, as computed according to MODE and COND, will
           be multiplied by a random sign (+1 or -1).  If RSIGN='F',
           they will not be.  RSIGN may only have the values 'T' or
           'F'.
           Not modified.
[in]UPPER
          UPPER is CHARACTER*1
           If UPPER='T', then the elements of A above the diagonal
           (and above the 2x2 diagonal blocks, if A has complex
           eigenvalues) will be set to random numbers out of DIST.
           If UPPER='F', they will not.  UPPER may only have the
           values 'T' or 'F'.
           Not modified.
[in]SIM
          SIM is CHARACTER*1
           If SIM='T', then A will be operated on by a "similarity
           transform", i.e., multiplied on the left by a matrix X and
           on the right by X inverse.  X = U S V, where U and V are
           random unitary matrices and S is a (diagonal) matrix of
           singular values specified by DS, MODES, and CONDS.  If
           SIM='F', then A will not be transformed.
           Not modified.
[in,out]DS
          DS is REAL array, dimension ( N )
           This array is used to specify the singular values of X,
           in the same way that D specifies the eigenvalues of A.
           If MODE=0, the DS contains the singular values, which
           may not be zero.
           Modified if MODE is nonzero.
[in]MODES
          MODES is INTEGER
[in]CONDS
          CONDS is REAL
           Same as MODE and COND, but for specifying the diagonal
           of S.  MODES=-6 and +6 are not allowed (since they would
           result in randomly ill-conditioned eigenvalues.)
[in]KL
          KL is INTEGER
           This specifies the lower bandwidth of the  matrix.  KL=1
           specifies upper Hessenberg form.  If KL is at least N-1,
           then A will have full lower bandwidth.  KL must be at
           least 1.
           Not modified.
[in]KU
          KU is INTEGER
           This specifies the upper bandwidth of the  matrix.  KU=1
           specifies lower Hessenberg form.  If KU is at least N-1,
           then A will have full upper bandwidth; if KU and KL
           are both at least N-1, then A will be dense.  Only one of
           KU and KL may be less than N-1.  KU must be at least 1.
           Not modified.
[in]ANORM
          ANORM is REAL
           If ANORM is not negative, then A will be scaled by a non-
           negative real number to make the maximum-element-norm of A
           to be ANORM.
           Not modified.
[out]A
          A is REAL array, dimension ( LDA, N )
           On exit A is the desired test matrix.
           Modified.
[in]LDA
          LDA is INTEGER
           LDA specifies the first dimension of A as declared in the
           calling program.  LDA must be at least N.
           Not modified.
[out]WORK
          WORK is REAL array, dimension ( 3*N )
           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 => N negative
            -2 => DIST illegal string
            -5 => MODE not in range -6 to 6
            -6 => COND less than 1.0, and MODE neither -6, 0 nor 6
            -8 => EI(1) is not ' ' or 'R', EI(j) is not 'R' or 'I', or
                  two adjacent elements of EI are 'I'.
            -9 => RSIGN is not 'T' or 'F'
           -10 => UPPER is not 'T' or 'F'
           -11 => SIM   is not 'T' or 'F'
           -12 => MODES=0 and DS has a zero singular value.
           -13 => MODES is not in the range -5 to 5.
           -14 => MODES is nonzero and CONDS is less than 1.
           -15 => KL is less than 1.
           -16 => KU is less than 1, or KL and KU are both less than
                  N-1.
           -19 => LDA is less than N.
            1  => Error return from SLATM1 (computing D)
            2  => Cannot scale to DMAX (max. eigenvalue is 0)
            3  => Error return from SLATM1 (computing DS)
            4  => Error return from SLARGE
            5  => Zero singular value from SLATM1.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 327 of file slatme.f.

332 *
333 * -- LAPACK computational routine --
334 * -- LAPACK is a software package provided by Univ. of Tennessee, --
335 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
336 *
337 * .. Scalar Arguments ..
338  CHARACTER DIST, RSIGN, SIM, UPPER
339  INTEGER INFO, KL, KU, LDA, MODE, MODES, N
340  REAL ANORM, COND, CONDS, DMAX
341 * ..
342 * .. Array Arguments ..
343  CHARACTER EI( * )
344  INTEGER ISEED( 4 )
345  REAL A( LDA, * ), D( * ), DS( * ), WORK( * )
346 * ..
347 *
348 * =====================================================================
349 *
350 * .. Parameters ..
351  REAL ZERO
352  parameter( zero = 0.0e0 )
353  REAL ONE
354  parameter( one = 1.0e0 )
355  REAL HALF
356  parameter( half = 1.0e0 / 2.0e0 )
357 * ..
358 * .. Local Scalars ..
359  LOGICAL BADEI, BADS, USEEI
360  INTEGER I, IC, ICOLS, IDIST, IINFO, IR, IROWS, IRSIGN,
361  $ ISIM, IUPPER, J, JC, JCR, JR
362  REAL ALPHA, TAU, TEMP, XNORMS
363 * ..
364 * .. Local Arrays ..
365  REAL TEMPA( 1 )
366 * ..
367 * .. External Functions ..
368  LOGICAL LSAME
369  REAL SLANGE, SLARAN
370  EXTERNAL lsame, slange, slaran
371 * ..
372 * .. External Subroutines ..
373  EXTERNAL scopy, sgemv, sger, slarfg, slarge, slarnv,
375 * ..
376 * .. Intrinsic Functions ..
377  INTRINSIC abs, max, mod
378 * ..
379 * .. Executable Statements ..
380 *
381 * 1) Decode and Test the input parameters.
382 * Initialize flags & seed.
383 *
384  info = 0
385 *
386 * Quick return if possible
387 *
388  IF( n.EQ.0 )
389  $ RETURN
390 *
391 * Decode DIST
392 *
393  IF( lsame( dist, 'U' ) ) THEN
394  idist = 1
395  ELSE IF( lsame( dist, 'S' ) ) THEN
396  idist = 2
397  ELSE IF( lsame( dist, 'N' ) ) THEN
398  idist = 3
399  ELSE
400  idist = -1
401  END IF
402 *
403 * Check EI
404 *
405  useei = .true.
406  badei = .false.
407  IF( lsame( ei( 1 ), ' ' ) .OR. mode.NE.0 ) THEN
408  useei = .false.
409  ELSE
410  IF( lsame( ei( 1 ), 'R' ) ) THEN
411  DO 10 j = 2, n
412  IF( lsame( ei( j ), 'I' ) ) THEN
413  IF( lsame( ei( j-1 ), 'I' ) )
414  $ badei = .true.
415  ELSE
416  IF( .NOT.lsame( ei( j ), 'R' ) )
417  $ badei = .true.
418  END IF
419  10 CONTINUE
420  ELSE
421  badei = .true.
422  END IF
423  END IF
424 *
425 * Decode RSIGN
426 *
427  IF( lsame( rsign, 'T' ) ) THEN
428  irsign = 1
429  ELSE IF( lsame( rsign, 'F' ) ) THEN
430  irsign = 0
431  ELSE
432  irsign = -1
433  END IF
434 *
435 * Decode UPPER
436 *
437  IF( lsame( upper, 'T' ) ) THEN
438  iupper = 1
439  ELSE IF( lsame( upper, 'F' ) ) THEN
440  iupper = 0
441  ELSE
442  iupper = -1
443  END IF
444 *
445 * Decode SIM
446 *
447  IF( lsame( sim, 'T' ) ) THEN
448  isim = 1
449  ELSE IF( lsame( sim, 'F' ) ) THEN
450  isim = 0
451  ELSE
452  isim = -1
453  END IF
454 *
455 * Check DS, if MODES=0 and ISIM=1
456 *
457  bads = .false.
458  IF( modes.EQ.0 .AND. isim.EQ.1 ) THEN
459  DO 20 j = 1, n
460  IF( ds( j ).EQ.zero )
461  $ bads = .true.
462  20 CONTINUE
463  END IF
464 *
465 * Set INFO if an error
466 *
467  IF( n.LT.0 ) THEN
468  info = -1
469  ELSE IF( idist.EQ.-1 ) THEN
470  info = -2
471  ELSE IF( abs( mode ).GT.6 ) THEN
472  info = -5
473  ELSE IF( ( mode.NE.0 .AND. abs( mode ).NE.6 ) .AND. cond.LT.one )
474  $ THEN
475  info = -6
476  ELSE IF( badei ) THEN
477  info = -8
478  ELSE IF( irsign.EQ.-1 ) THEN
479  info = -9
480  ELSE IF( iupper.EQ.-1 ) THEN
481  info = -10
482  ELSE IF( isim.EQ.-1 ) THEN
483  info = -11
484  ELSE IF( bads ) THEN
485  info = -12
486  ELSE IF( isim.EQ.1 .AND. abs( modes ).GT.5 ) THEN
487  info = -13
488  ELSE IF( isim.EQ.1 .AND. modes.NE.0 .AND. conds.LT.one ) THEN
489  info = -14
490  ELSE IF( kl.LT.1 ) THEN
491  info = -15
492  ELSE IF( ku.LT.1 .OR. ( ku.LT.n-1 .AND. kl.LT.n-1 ) ) THEN
493  info = -16
494  ELSE IF( lda.LT.max( 1, n ) ) THEN
495  info = -19
496  END IF
497 *
498  IF( info.NE.0 ) THEN
499  CALL xerbla( 'SLATME', -info )
500  RETURN
501  END IF
502 *
503 * Initialize random number generator
504 *
505  DO 30 i = 1, 4
506  iseed( i ) = mod( abs( iseed( i ) ), 4096 )
507  30 CONTINUE
508 *
509  IF( mod( iseed( 4 ), 2 ).NE.1 )
510  $ iseed( 4 ) = iseed( 4 ) + 1
511 *
512 * 2) Set up diagonal of A
513 *
514 * Compute D according to COND and MODE
515 *
516  CALL slatm1( mode, cond, irsign, idist, iseed, d, n, iinfo )
517  IF( iinfo.NE.0 ) THEN
518  info = 1
519  RETURN
520  END IF
521  IF( mode.NE.0 .AND. abs( mode ).NE.6 ) THEN
522 *
523 * Scale by DMAX
524 *
525  temp = abs( d( 1 ) )
526  DO 40 i = 2, n
527  temp = max( temp, abs( d( i ) ) )
528  40 CONTINUE
529 *
530  IF( temp.GT.zero ) THEN
531  alpha = dmax / temp
532  ELSE IF( dmax.NE.zero ) THEN
533  info = 2
534  RETURN
535  ELSE
536  alpha = zero
537  END IF
538 *
539  CALL sscal( n, alpha, d, 1 )
540 *
541  END IF
542 *
543  CALL slaset( 'Full', n, n, zero, zero, a, lda )
544  CALL scopy( n, d, 1, a, lda+1 )
545 *
546 * Set up complex conjugate pairs
547 *
548  IF( mode.EQ.0 ) THEN
549  IF( useei ) THEN
550  DO 50 j = 2, n
551  IF( lsame( ei( j ), 'I' ) ) THEN
552  a( j-1, j ) = a( j, j )
553  a( j, j-1 ) = -a( j, j )
554  a( j, j ) = a( j-1, j-1 )
555  END IF
556  50 CONTINUE
557  END IF
558 *
559  ELSE IF( abs( mode ).EQ.5 ) THEN
560 *
561  DO 60 j = 2, n, 2
562  IF( slaran( iseed ).GT.half ) THEN
563  a( j-1, j ) = a( j, j )
564  a( j, j-1 ) = -a( j, j )
565  a( j, j ) = a( j-1, j-1 )
566  END IF
567  60 CONTINUE
568  END IF
569 *
570 * 3) If UPPER='T', set upper triangle of A to random numbers.
571 * (but don't modify the corners of 2x2 blocks.)
572 *
573  IF( iupper.NE.0 ) THEN
574  DO 70 jc = 2, n
575  IF( a( jc-1, jc ).NE.zero ) THEN
576  jr = jc - 2
577  ELSE
578  jr = jc - 1
579  END IF
580  CALL slarnv( idist, iseed, jr, a( 1, jc ) )
581  70 CONTINUE
582  END IF
583 *
584 * 4) If SIM='T', apply similarity transformation.
585 *
586 * -1
587 * Transform is X A X , where X = U S V, thus
588 *
589 * it is U S V A V' (1/S) U'
590 *
591  IF( isim.NE.0 ) THEN
592 *
593 * Compute S (singular values of the eigenvector matrix)
594 * according to CONDS and MODES
595 *
596  CALL slatm1( modes, conds, 0, 0, iseed, ds, n, iinfo )
597  IF( iinfo.NE.0 ) THEN
598  info = 3
599  RETURN
600  END IF
601 *
602 * Multiply by V and V'
603 *
604  CALL slarge( n, a, lda, iseed, work, iinfo )
605  IF( iinfo.NE.0 ) THEN
606  info = 4
607  RETURN
608  END IF
609 *
610 * Multiply by S and (1/S)
611 *
612  DO 80 j = 1, n
613  CALL sscal( n, ds( j ), a( j, 1 ), lda )
614  IF( ds( j ).NE.zero ) THEN
615  CALL sscal( n, one / ds( j ), a( 1, j ), 1 )
616  ELSE
617  info = 5
618  RETURN
619  END IF
620  80 CONTINUE
621 *
622 * Multiply by U and U'
623 *
624  CALL slarge( n, a, lda, iseed, work, iinfo )
625  IF( iinfo.NE.0 ) THEN
626  info = 4
627  RETURN
628  END IF
629  END IF
630 *
631 * 5) Reduce the bandwidth.
632 *
633  IF( kl.LT.n-1 ) THEN
634 *
635 * Reduce bandwidth -- kill column
636 *
637  DO 90 jcr = kl + 1, n - 1
638  ic = jcr - kl
639  irows = n + 1 - jcr
640  icols = n + kl - jcr
641 *
642  CALL scopy( irows, a( jcr, ic ), 1, work, 1 )
643  xnorms = work( 1 )
644  CALL slarfg( irows, xnorms, work( 2 ), 1, tau )
645  work( 1 ) = one
646 *
647  CALL sgemv( 'T', irows, icols, one, a( jcr, ic+1 ), lda,
648  $ work, 1, zero, work( irows+1 ), 1 )
649  CALL sger( irows, icols, -tau, work, 1, work( irows+1 ), 1,
650  $ a( jcr, ic+1 ), lda )
651 *
652  CALL sgemv( 'N', n, irows, one, a( 1, jcr ), lda, work, 1,
653  $ zero, work( irows+1 ), 1 )
654  CALL sger( n, irows, -tau, work( irows+1 ), 1, work, 1,
655  $ a( 1, jcr ), lda )
656 *
657  a( jcr, ic ) = xnorms
658  CALL slaset( 'Full', irows-1, 1, zero, zero, a( jcr+1, ic ),
659  $ lda )
660  90 CONTINUE
661  ELSE IF( ku.LT.n-1 ) THEN
662 *
663 * Reduce upper bandwidth -- kill a row at a time.
664 *
665  DO 100 jcr = ku + 1, n - 1
666  ir = jcr - ku
667  irows = n + ku - jcr
668  icols = n + 1 - jcr
669 *
670  CALL scopy( icols, a( ir, jcr ), lda, work, 1 )
671  xnorms = work( 1 )
672  CALL slarfg( icols, xnorms, work( 2 ), 1, tau )
673  work( 1 ) = one
674 *
675  CALL sgemv( 'N', irows, icols, one, a( ir+1, jcr ), lda,
676  $ work, 1, zero, work( icols+1 ), 1 )
677  CALL sger( irows, icols, -tau, work( icols+1 ), 1, work, 1,
678  $ a( ir+1, jcr ), lda )
679 *
680  CALL sgemv( 'C', icols, n, one, a( jcr, 1 ), lda, work, 1,
681  $ zero, work( icols+1 ), 1 )
682  CALL sger( icols, n, -tau, work, 1, work( icols+1 ), 1,
683  $ a( jcr, 1 ), lda )
684 *
685  a( ir, jcr ) = xnorms
686  CALL slaset( 'Full', 1, icols-1, zero, zero, a( ir, jcr+1 ),
687  $ lda )
688  100 CONTINUE
689  END IF
690 *
691 * Scale the matrix to have norm ANORM
692 *
693  IF( anorm.GE.zero ) THEN
694  temp = slange( 'M', n, n, a, lda, tempa )
695  IF( temp.GT.zero ) THEN
696  alpha = anorm / temp
697  DO 110 j = 1, n
698  CALL sscal( n, alpha, a( 1, j ), 1 )
699  110 CONTINUE
700  END IF
701  END IF
702 *
703  RETURN
704 *
705 * End of SLATME
706 *
subroutine slarnv(IDIST, ISEED, N, X)
SLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition: slarnv.f:97
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 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
real function slaran(ISEED)
SLARAN
Definition: slaran.f:67
subroutine slarge(N, A, LDA, ISEED, WORK, INFO)
SLARGE
Definition: slarge.f:87
real function slange(NORM, M, N, A, LDA, WORK)
SLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: slange.f:114
subroutine slarfg(N, ALPHA, X, INCX, TAU)
SLARFG generates an elementary reflector (Householder matrix).
Definition: slarfg.f:106
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:82
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:79
subroutine sger(M, N, ALPHA, X, INCX, Y, INCY, A, LDA)
SGER
Definition: sger.f:130
subroutine sgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
SGEMV
Definition: sgemv.f:156
Here is the call graph for this function:
Here is the caller graph for this function: