LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ clatmr()

subroutine clatmr ( integer m,
integer n,
character dist,
integer, dimension( 4 ) iseed,
character sym,
complex, dimension( * ) d,
integer mode,
real cond,
complex dmax,
character rsign,
character grade,
complex, dimension( * ) dl,
integer model,
real condl,
complex, dimension( * ) dr,
integer moder,
real condr,
character pivtng,
integer, dimension( * ) ipivot,
integer kl,
integer ku,
real sparse,
real anorm,
character pack,
complex, dimension( lda, * ) a,
integer lda,
integer, dimension( * ) iwork,
integer info )

CLATMR

Purpose:
!>
!>    CLATMR generates random matrices of various types for testing
!>    LAPACK programs.
!>
!>    CLATMR operates by applying the following sequence of
!>    operations:
!>
!>      Generate a matrix A with random entries of distribution DIST
!>         which is symmetric if SYM='S', Hermitian if SYM='H', and
!>         nonsymmetric if SYM='N'.
!>
!>      Set the diagonal to D, where D may be input or
!>         computed according to MODE, COND, DMAX and RSIGN
!>         as described below.
!>
!>      Grade the matrix, if desired, from the left and/or right
!>         as specified by GRADE. The inputs DL, MODEL, CONDL, DR,
!>         MODER and CONDR also determine the grading as described
!>         below.
!>
!>      Permute, if desired, the rows and/or columns as specified by
!>         PIVTNG and IPIVOT.
!>
!>      Set random entries to zero, if desired, to get a random sparse
!>         matrix as specified by SPARSE.
!>
!>      Make A a band matrix, if desired, by zeroing out the matrix
!>         outside a band of lower bandwidth KL and upper bandwidth KU.
!>
!>      Scale A, if desired, to have maximum entry ANORM.
!>
!>      Pack the matrix if desired. Options specified by PACK are:
!>         no packing
!>         zero out upper half (if symmetric or Hermitian)
!>         zero out lower half (if symmetric or Hermitian)
!>         store the upper half columnwise (if symmetric or Hermitian
!>             or square upper triangular)
!>         store the lower half columnwise (if symmetric or Hermitian
!>             or square lower triangular)
!>             same as upper half rowwise if symmetric
!>             same as conjugate upper half rowwise if Hermitian
!>         store the lower triangle in banded format
!>             (if symmetric or Hermitian)
!>         store the upper triangle in banded format
!>             (if symmetric or Hermitian)
!>         store the entire matrix in banded format
!>
!>    Note: If two calls to CLATMR differ only in the PACK parameter,
!>          they will generate mathematically equivalent matrices.
!>
!>          If two calls to CLATMR both have full bandwidth (KL = M-1
!>          and KU = N-1), and differ only in the PIVTNG and PACK
!>          parameters, then the matrices generated will differ only
!>          in the order of the rows and/or columns, and otherwise
!>          contain the same data. This consistency cannot be and
!>          is not maintained with less than full bandwidth.
!> 
Parameters
[in]M
!>          M is INTEGER
!>           Number of rows of A. Not modified.
!> 
[in]N
!>          N is INTEGER
!>           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 a random matrix .
!>           'U' => real and imaginary parts are independent
!>                  UNIFORM( 0, 1 )  ( 'U' for uniform )
!>           'S' => real and imaginary parts are independent
!>                  UNIFORM( -1, 1 ) ( 'S' for symmetric )
!>           'N' => real and imaginary parts are independent
!>                  NORMAL( 0, 1 )   ( 'N' for normal )
!>           'D' => uniform on interior of unit disk ( 'D' for disk )
!>           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 CLATMR
!>           to continue the same random number sequence.
!>           Changed on exit.
!> 
[in]SYM
!>          SYM is CHARACTER*1
!>           If SYM='S', generated matrix is symmetric.
!>           If SYM='H', generated matrix is Hermitian.
!>           If SYM='N', generated matrix is nonsymmetric.
!>           Not modified.
!> 
[in,out]D
!>          D is COMPLEX array, dimension (min(M,N))
!>           On entry this array specifies the diagonal entries
!>           of the diagonal of A.  D may either be specified
!>           on entry, or set according to MODE and COND as described
!>           below. If the matrix is Hermitian, the real part of D
!>           will be taken. May be changed on exit if MODE is nonzero.
!> 
[in]MODE
!>          MODE is INTEGER
!>           On entry describes how D is to be used:
!>           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,
!>           Not modified.
!> 
[in]COND
!>          COND is REAL
!>           On entry, used as described under MODE above.
!>           If used, it must be >= 1. Not modified.
!> 
[in]DMAX
!>          DMAX is COMPLEX
!>           If MODE neither -6, 0 nor 6, the diagonal is scaled by
!>           DMAX / max(abs(D(i))), so that maximum absolute entry
!>           of diagonal is abs(DMAX). If DMAX is complex (or zero),
!>           diagonal will be scaled by a complex number (or zero).
!> 
[in]RSIGN
!>          RSIGN is CHARACTER*1
!>           If MODE neither -6, 0 nor 6, specifies sign of diagonal
!>           as follows:
!>           'T' => diagonal entries are multiplied by a random complex
!>                  number uniformly distributed with absolute value 1
!>           'F' => diagonal unchanged
!>           Not modified.
!> 
[in]GRADE
!>          GRADE is CHARACTER*1
!>           Specifies grading of matrix as follows:
!>           'N'  => no grading
!>           'L'  => matrix premultiplied by diag( DL )
!>                   (only if matrix nonsymmetric)
!>           'R'  => matrix postmultiplied by diag( DR )
!>                   (only if matrix nonsymmetric)
!>           'B'  => matrix premultiplied by diag( DL ) and
!>                         postmultiplied by diag( DR )
!>                   (only if matrix nonsymmetric)
!>           'H'  => matrix premultiplied by diag( DL ) and
!>                         postmultiplied by diag( CONJG(DL) )
!>                   (only if matrix Hermitian or nonsymmetric)
!>           'S'  => matrix premultiplied by diag( DL ) and
!>                         postmultiplied by diag( DL )
!>                   (only if matrix symmetric or nonsymmetric)
!>           'E'  => matrix premultiplied by diag( DL ) and
!>                         postmultiplied by inv( diag( DL ) )
!>                         ( 'S' for similarity )
!>                   (only if matrix nonsymmetric)
!>                   Note: if GRADE='S', then M must equal N.
!>           Not modified.
!> 
[in,out]DL
!>          DL is COMPLEX array, dimension (M)
!>           If MODEL=0, then on entry this array specifies the diagonal
!>           entries of a diagonal matrix used as described under GRADE
!>           above. If MODEL is not zero, then DL will be set according
!>           to MODEL and CONDL, analogous to the way D is set according
!>           to MODE and COND (except there is no DMAX parameter for DL).
!>           If GRADE='E', then DL cannot have zero entries.
!>           Not referenced if GRADE = 'N' or 'R'. Changed on exit.
!> 
[in]MODEL
!>          MODEL is INTEGER
!>           This specifies how the diagonal array DL is to be computed,
!>           just as MODE specifies how D is to be computed.
!>           Not modified.
!> 
[in]CONDL
!>          CONDL is REAL
!>           When MODEL is not zero, this specifies the condition number
!>           of the computed DL.  Not modified.
!> 
[in,out]DR
!>          DR is COMPLEX array, dimension (N)
!>           If MODER=0, then on entry this array specifies the diagonal
!>           entries of a diagonal matrix used as described under GRADE
!>           above. If MODER is not zero, then DR will be set according
!>           to MODER and CONDR, analogous to the way D is set according
!>           to MODE and COND (except there is no DMAX parameter for DR).
!>           Not referenced if GRADE = 'N', 'L', 'H' or 'S'.
!>           Changed on exit.
!> 
[in]MODER
!>          MODER is INTEGER
!>           This specifies how the diagonal array DR is to be computed,
!>           just as MODE specifies how D is to be computed.
!>           Not modified.
!> 
[in]CONDR
!>          CONDR is REAL
!>           When MODER is not zero, this specifies the condition number
!>           of the computed DR.  Not modified.
!> 
[in]PIVTNG
!>          PIVTNG is CHARACTER*1
!>           On entry specifies pivoting permutations as follows:
!>           'N' or ' ' => none.
!>           'L' => left or row pivoting (matrix must be nonsymmetric).
!>           'R' => right or column pivoting (matrix must be
!>                  nonsymmetric).
!>           'B' or 'F' => both or full pivoting, i.e., on both sides.
!>                         In this case, M must equal N
!>
!>           If two calls to CLATMR both have full bandwidth (KL = M-1
!>           and KU = N-1), and differ only in the PIVTNG and PACK
!>           parameters, then the matrices generated will differ only
!>           in the order of the rows and/or columns, and otherwise
!>           contain the same data. This consistency cannot be
!>           maintained with less than full bandwidth.
!> 
[in]IPIVOT
!>          IPIVOT is INTEGER array, dimension (N or M)
!>           This array specifies the permutation used.  After the
!>           basic matrix is generated, the rows, columns, or both
!>           are permuted.   If, say, row pivoting is selected, CLATMR
!>           starts with the *last* row and interchanges the M-th and
!>           IPIVOT(M)-th rows, then moves to the next-to-last row,
!>           interchanging the (M-1)-th and the IPIVOT(M-1)-th rows,
!>           and so on.  In terms of , the permutation is
!>           (1 IPIVOT(1)) (2 IPIVOT(2)) ... (M IPIVOT(M))
!>           where the rightmost cycle is applied first.  This is the
!>           *inverse* of the effect of pivoting in LINPACK.  The idea
!>           is that factoring (with pivoting) an identity matrix
!>           which has been inverse-pivoted in this way should
!>           result in a pivot vector identical to IPIVOT.
!>           Not referenced if PIVTNG = 'N'. Not modified.
!> 
[in]KL
!>          KL is INTEGER
!>           On entry specifies the lower bandwidth of the  matrix. For
!>           example, KL=0 implies upper triangular, KL=1 implies upper
!>           Hessenberg, and KL at least M-1 implies the matrix is not
!>           banded. Must equal KU if matrix is symmetric or Hermitian.
!>           Not modified.
!> 
[in]KU
!>          KU is INTEGER
!>           On entry specifies the upper bandwidth of the  matrix. For
!>           example, KU=0 implies lower triangular, KU=1 implies lower
!>           Hessenberg, and KU at least N-1 implies the matrix is not
!>           banded. Must equal KL if matrix is symmetric or Hermitian.
!>           Not modified.
!> 
[in]SPARSE
!>          SPARSE is REAL
!>           On entry specifies the sparsity of the matrix if a sparse
!>           matrix is to be generated. SPARSE should lie between
!>           0 and 1. To generate a sparse matrix, for each matrix entry
!>           a uniform ( 0, 1 ) random number x is generated and
!>           compared to SPARSE; if x is larger the matrix entry
!>           is unchanged and if x is smaller the entry is set
!>           to zero. Thus on the average a fraction SPARSE of the
!>           entries will be set to zero.
!>           Not modified.
!> 
[in]ANORM
!>          ANORM is REAL
!>           On entry specifies maximum entry of output matrix
!>           (output matrix will by multiplied by a constant so that
!>           its largest absolute entry equal ANORM)
!>           if ANORM is nonnegative. If ANORM is negative no scaling
!>           is done. Not modified.
!> 
[in]PACK
!>          PACK is CHARACTER*1
!>           On entry specifies packing of matrix as follows:
!>           'N' => no packing
!>           'U' => zero out all subdiagonal entries
!>                  (if symmetric or Hermitian)
!>           'L' => zero out all superdiagonal entries
!>                  (if symmetric or Hermitian)
!>           'C' => store the upper triangle columnwise
!>                  (only if matrix symmetric or Hermitian or
!>                   square upper triangular)
!>           'R' => store the lower triangle columnwise
!>                  (only if matrix symmetric or Hermitian or
!>                   square lower triangular)
!>                  (same as upper half rowwise if symmetric)
!>                  (same as conjugate upper half rowwise if Hermitian)
!>           'B' => store the lower triangle in band storage scheme
!>                  (only if matrix symmetric or Hermitian)
!>           'Q' => store the upper triangle in band storage scheme
!>                  (only if matrix symmetric or Hermitian)
!>           '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, HB or TB     - use 'B' or 'Q'
!>           PP, HP or TP     - use 'C' or 'R'
!>
!>           If two calls to CLATMR differ only in the PACK parameter,
!>           they will generate mathematically equivalent matrices.
!>           Not modified.
!> 
[in,out]A
!>          A is COMPLEX array, dimension (LDA,N)
!>           On exit A is the desired test matrix. Only those
!>           entries of A which are significant on output
!>           will be referenced (even if A is in packed or band
!>           storage format). The 'unoccupied corners' of A in
!>           band format will be zeroed out.
!> 
[in]LDA
!>          LDA is INTEGER
!>           on entry LDA specifies the first dimension of A as
!>           declared in the calling program.
!>           If PACK='N', 'U' or 'L', LDA must be at least max ( 1, M ).
!>           If PACK='C' or 'R', LDA must be at least 1.
!>           If PACK='B', or 'Q', LDA must be MIN ( KU+1, N )
!>           If PACK='Z', LDA must be at least KUU+KLL+1, where
!>           KUU = MIN ( KU, N-1 ) and KLL = MIN ( KL, M-1 )
!>           Not modified.
!> 
[out]IWORK
!>          IWORK is INTEGER array, dimension (N or M)
!>           Workspace. Not referenced if PIVTNG = 'N'. Changed on exit.
!> 
[out]INFO
!>          INFO is INTEGER
!>           Error parameter on exit:
!>             0 => normal return
!>            -1 => M negative or unequal to N and SYM='S' or 'H'
!>            -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 => MODE neither -6, 0 nor 6 and RSIGN illegal string
!>           -11 => GRADE illegal string, or GRADE='E' and
!>                  M not equal to N, or GRADE='L', 'R', 'B', 'S' or 'E'
!>                  and SYM = 'H', or GRADE='L', 'R', 'B', 'H' or 'E'
!>                  and SYM = 'S'
!>           -12 => GRADE = 'E' and DL contains zero
!>           -13 => MODEL not in range -6 to 6 and GRADE= 'L', 'B', 'H',
!>                  'S' or 'E'
!>           -14 => CONDL less than 1.0, GRADE='L', 'B', 'H', 'S' or 'E',
!>                  and MODEL neither -6, 0 nor 6
!>           -16 => MODER not in range -6 to 6 and GRADE= 'R' or 'B'
!>           -17 => CONDR less than 1.0, GRADE='R' or 'B', and
!>                  MODER neither -6, 0 nor 6
!>           -18 => PIVTNG illegal string, or PIVTNG='B' or 'F' and
!>                  M not equal to N, or PIVTNG='L' or 'R' and SYM='S'
!>                  or 'H'
!>           -19 => IPIVOT contains out of range number and
!>                  PIVTNG not equal to 'N'
!>           -20 => KL negative
!>           -21 => KU negative, or SYM='S' or 'H' and KU not equal to KL
!>           -22 => SPARSE not in range 0. to 1.
!>           -24 => PACK illegal string, or PACK='U', 'L', 'B' or 'Q'
!>                  and SYM='N', or PACK='C' and SYM='N' and either KL
!>                  not equal to 0 or N not equal to M, or PACK='R' and
!>                  SYM='N', and either KU not equal to 0 or N not equal
!>                  to M
!>           -26 => LDA too small
!>             1 => Error return from CLATM1 (computing D)
!>             2 => Cannot scale diagonal to DMAX (max. entry is 0)
!>             3 => Error return from CLATM1 (computing DL)
!>             4 => Error return from CLATM1 (computing DR)
!>             5 => ANORM is positive, but matrix constructed prior to
!>                  attempting to scale it to have norm ANORM, is zero
!> 
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 486 of file clatmr.f.

490*
491* -- LAPACK computational routine --
492* -- LAPACK is a software package provided by Univ. of Tennessee, --
493* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
494*
495* .. Scalar Arguments ..
496 CHARACTER DIST, GRADE, PACK, PIVTNG, RSIGN, SYM
497 INTEGER INFO, KL, KU, LDA, M, MODE, MODEL, MODER, N
498 REAL ANORM, COND, CONDL, CONDR, SPARSE
499 COMPLEX DMAX
500* ..
501* .. Array Arguments ..
502 INTEGER IPIVOT( * ), ISEED( 4 ), IWORK( * )
503 COMPLEX A( LDA, * ), D( * ), DL( * ), DR( * )
504* ..
505*
506* =====================================================================
507*
508* .. Parameters ..
509 REAL ZERO
510 parameter( zero = 0.0e0 )
511 REAL ONE
512 parameter( one = 1.0e0 )
513 COMPLEX CONE
514 parameter( cone = ( 1.0e0, 0.0e0 ) )
515 COMPLEX CZERO
516 parameter( czero = ( 0.0e0, 0.0e0 ) )
517* ..
518* .. Local Scalars ..
519 LOGICAL BADPVT, DZERO, FULBND
520 INTEGER I, IDIST, IGRADE, IISUB, IPACK, IPVTNG, IRSIGN,
521 $ ISUB, ISYM, J, JJSUB, JSUB, K, KLL, KUU, MNMIN,
522 $ MNSUB, MXSUB, NPVTS
523 REAL ONORM, TEMP
524 COMPLEX CALPHA, CTEMP
525* ..
526* .. Local Arrays ..
527 REAL TEMPA( 1 )
528* ..
529* .. External Functions ..
530 LOGICAL LSAME
531 REAL CLANGB, CLANGE, CLANSB, CLANSP, CLANSY
532 COMPLEX CLATM2, CLATM3
533 EXTERNAL lsame, clangb, clange,
534 $ clansb, clansp, clansy,
535 $ clatm2, clatm3
536* ..
537* .. External Subroutines ..
538 EXTERNAL clatm1, csscal, xerbla
539* ..
540* .. Intrinsic Functions ..
541 INTRINSIC abs, conjg, max, min, mod, real
542* ..
543* .. Executable Statements ..
544*
545* 1) Decode and Test the input parameters.
546* Initialize flags & seed.
547*
548 info = 0
549*
550* Quick return if possible
551*
552 IF( m.EQ.0 .OR. n.EQ.0 )
553 $ RETURN
554*
555* Decode DIST
556*
557 IF( lsame( dist, 'U' ) ) THEN
558 idist = 1
559 ELSE IF( lsame( dist, 'S' ) ) THEN
560 idist = 2
561 ELSE IF( lsame( dist, 'N' ) ) THEN
562 idist = 3
563 ELSE IF( lsame( dist, 'D' ) ) THEN
564 idist = 4
565 ELSE
566 idist = -1
567 END IF
568*
569* Decode SYM
570*
571 IF( lsame( sym, 'H' ) ) THEN
572 isym = 0
573 ELSE IF( lsame( sym, 'N' ) ) THEN
574 isym = 1
575 ELSE IF( lsame( sym, 'S' ) ) THEN
576 isym = 2
577 ELSE
578 isym = -1
579 END IF
580*
581* Decode RSIGN
582*
583 IF( lsame( rsign, 'F' ) ) THEN
584 irsign = 0
585 ELSE IF( lsame( rsign, 'T' ) ) THEN
586 irsign = 1
587 ELSE
588 irsign = -1
589 END IF
590*
591* Decode PIVTNG
592*
593 IF( lsame( pivtng, 'N' ) ) THEN
594 ipvtng = 0
595 ELSE IF( lsame( pivtng, ' ' ) ) THEN
596 ipvtng = 0
597 ELSE IF( lsame( pivtng, 'L' ) ) THEN
598 ipvtng = 1
599 npvts = m
600 ELSE IF( lsame( pivtng, 'R' ) ) THEN
601 ipvtng = 2
602 npvts = n
603 ELSE IF( lsame( pivtng, 'B' ) ) THEN
604 ipvtng = 3
605 npvts = min( n, m )
606 ELSE IF( lsame( pivtng, 'F' ) ) THEN
607 ipvtng = 3
608 npvts = min( n, m )
609 ELSE
610 ipvtng = -1
611 END IF
612*
613* Decode GRADE
614*
615 IF( lsame( grade, 'N' ) ) THEN
616 igrade = 0
617 ELSE IF( lsame( grade, 'L' ) ) THEN
618 igrade = 1
619 ELSE IF( lsame( grade, 'R' ) ) THEN
620 igrade = 2
621 ELSE IF( lsame( grade, 'B' ) ) THEN
622 igrade = 3
623 ELSE IF( lsame( grade, 'E' ) ) THEN
624 igrade = 4
625 ELSE IF( lsame( grade, 'H' ) ) THEN
626 igrade = 5
627 ELSE IF( lsame( grade, 'S' ) ) THEN
628 igrade = 6
629 ELSE
630 igrade = -1
631 END IF
632*
633* Decode PACK
634*
635 IF( lsame( pack, 'N' ) ) THEN
636 ipack = 0
637 ELSE IF( lsame( pack, 'U' ) ) THEN
638 ipack = 1
639 ELSE IF( lsame( pack, 'L' ) ) THEN
640 ipack = 2
641 ELSE IF( lsame( pack, 'C' ) ) THEN
642 ipack = 3
643 ELSE IF( lsame( pack, 'R' ) ) THEN
644 ipack = 4
645 ELSE IF( lsame( pack, 'B' ) ) THEN
646 ipack = 5
647 ELSE IF( lsame( pack, 'Q' ) ) THEN
648 ipack = 6
649 ELSE IF( lsame( pack, 'Z' ) ) THEN
650 ipack = 7
651 ELSE
652 ipack = -1
653 END IF
654*
655* Set certain internal parameters
656*
657 mnmin = min( m, n )
658 kll = min( kl, m-1 )
659 kuu = min( ku, n-1 )
660*
661* If inv(DL) is used, check to see if DL has a zero entry.
662*
663 dzero = .false.
664 IF( igrade.EQ.4 .AND. model.EQ.0 ) THEN
665 DO 10 i = 1, m
666 IF( dl( i ).EQ.czero )
667 $ dzero = .true.
668 10 CONTINUE
669 END IF
670*
671* Check values in IPIVOT
672*
673 badpvt = .false.
674 IF( ipvtng.GT.0 ) THEN
675 DO 20 j = 1, npvts
676 IF( ipivot( j ).LE.0 .OR. ipivot( j ).GT.npvts )
677 $ badpvt = .true.
678 20 CONTINUE
679 END IF
680*
681* Set INFO if an error
682*
683 IF( m.LT.0 ) THEN
684 info = -1
685 ELSE IF( m.NE.n .AND. ( isym.EQ.0 .OR. isym.EQ.2 ) ) THEN
686 info = -1
687 ELSE IF( n.LT.0 ) THEN
688 info = -2
689 ELSE IF( idist.EQ.-1 ) THEN
690 info = -3
691 ELSE IF( isym.EQ.-1 ) THEN
692 info = -5
693 ELSE IF( mode.LT.-6 .OR. mode.GT.6 ) THEN
694 info = -7
695 ELSE IF( ( mode.NE.-6 .AND. mode.NE.0 .AND. mode.NE.6 ) .AND.
696 $ cond.LT.one ) THEN
697 info = -8
698 ELSE IF( ( mode.NE.-6 .AND. mode.NE.0 .AND. mode.NE.6 ) .AND.
699 $ irsign.EQ.-1 ) THEN
700 info = -10
701 ELSE IF( igrade.EQ.-1 .OR. ( igrade.EQ.4 .AND. m.NE.n ) .OR.
702 $ ( ( igrade.EQ.1 .OR. igrade.EQ.2 .OR. igrade.EQ.3 .OR.
703 $ igrade.EQ.4 .OR. igrade.EQ.6 ) .AND. isym.EQ.0 ) .OR.
704 $ ( ( igrade.EQ.1 .OR. igrade.EQ.2 .OR. igrade.EQ.3 .OR.
705 $ igrade.EQ.4 .OR. igrade.EQ.5 ) .AND. isym.EQ.2 ) ) THEN
706 info = -11
707 ELSE IF( igrade.EQ.4 .AND. dzero ) THEN
708 info = -12
709 ELSE IF( ( igrade.EQ.1 .OR. igrade.EQ.3 .OR. igrade.EQ.4 .OR.
710 $ igrade.EQ.5 .OR. igrade.EQ.6 ) .AND.
711 $ ( model.LT.-6 .OR. model.GT.6 ) ) THEN
712 info = -13
713 ELSE IF( ( igrade.EQ.1 .OR. igrade.EQ.3 .OR. igrade.EQ.4 .OR.
714 $ igrade.EQ.5 .OR. igrade.EQ.6 ) .AND.
715 $ ( model.NE.-6 .AND. model.NE.0 .AND. model.NE.6 ) .AND.
716 $ condl.LT.one ) THEN
717 info = -14
718 ELSE IF( ( igrade.EQ.2 .OR. igrade.EQ.3 ) .AND.
719 $ ( moder.LT.-6 .OR. moder.GT.6 ) ) THEN
720 info = -16
721 ELSE IF( ( igrade.EQ.2 .OR. igrade.EQ.3 ) .AND.
722 $ ( moder.NE.-6 .AND. moder.NE.0 .AND. moder.NE.6 ) .AND.
723 $ condr.LT.one ) THEN
724 info = -17
725 ELSE IF( ipvtng.EQ.-1 .OR. ( ipvtng.EQ.3 .AND. m.NE.n ) .OR.
726 $ ( ( ipvtng.EQ.1 .OR. ipvtng.EQ.2 ) .AND. ( isym.EQ.0 .OR.
727 $ isym.EQ.2 ) ) ) THEN
728 info = -18
729 ELSE IF( ipvtng.NE.0 .AND. badpvt ) THEN
730 info = -19
731 ELSE IF( kl.LT.0 ) THEN
732 info = -20
733 ELSE IF( ku.LT.0 .OR. ( ( isym.EQ.0 .OR. isym.EQ.2 ) .AND. kl.NE.
734 $ ku ) ) THEN
735 info = -21
736 ELSE IF( sparse.LT.zero .OR. sparse.GT.one ) THEN
737 info = -22
738 ELSE IF( ipack.EQ.-1 .OR. ( ( ipack.EQ.1 .OR. ipack.EQ.2 .OR.
739 $ ipack.EQ.5 .OR. ipack.EQ.6 ) .AND. isym.EQ.1 ) .OR.
740 $ ( ipack.EQ.3 .AND. isym.EQ.1 .AND. ( kl.NE.0 .OR. m.NE.
741 $ n ) ) .OR. ( ipack.EQ.4 .AND. isym.EQ.1 .AND. ( ku.NE.
742 $ 0 .OR. m.NE.n ) ) ) THEN
743 info = -24
744 ELSE IF( ( ( ipack.EQ.0 .OR. ipack.EQ.1 .OR. ipack.EQ.2 ) .AND.
745 $ lda.LT.max( 1, m ) ) .OR. ( ( ipack.EQ.3 .OR. ipack.EQ.
746 $ 4 ) .AND. lda.LT.1 ) .OR. ( ( ipack.EQ.5 .OR. ipack.EQ.
747 $ 6 ) .AND. lda.LT.kuu+1 ) .OR.
748 $ ( ipack.EQ.7 .AND. lda.LT.kll+kuu+1 ) ) THEN
749 info = -26
750 END IF
751*
752 IF( info.NE.0 ) THEN
753 CALL xerbla( 'CLATMR', -info )
754 RETURN
755 END IF
756*
757* Decide if we can pivot consistently
758*
759 fulbnd = .false.
760 IF( kuu.EQ.n-1 .AND. kll.EQ.m-1 )
761 $ fulbnd = .true.
762*
763* Initialize random number generator
764*
765 DO 30 i = 1, 4
766 iseed( i ) = mod( abs( iseed( i ) ), 4096 )
767 30 CONTINUE
768*
769 iseed( 4 ) = 2*( iseed( 4 ) / 2 ) + 1
770*
771* 2) Set up D, DL, and DR, if indicated.
772*
773* Compute D according to COND and MODE
774*
775 CALL clatm1( mode, cond, irsign, idist, iseed, d, mnmin, info )
776 IF( info.NE.0 ) THEN
777 info = 1
778 RETURN
779 END IF
780 IF( mode.NE.0 .AND. mode.NE.-6 .AND. mode.NE.6 ) THEN
781*
782* Scale by DMAX
783*
784 temp = abs( d( 1 ) )
785 DO 40 i = 2, mnmin
786 temp = max( temp, abs( d( i ) ) )
787 40 CONTINUE
788 IF( temp.EQ.zero .AND. dmax.NE.czero ) THEN
789 info = 2
790 RETURN
791 END IF
792 IF( temp.NE.zero ) THEN
793 calpha = dmax / temp
794 ELSE
795 calpha = cone
796 END IF
797 DO 50 i = 1, mnmin
798 d( i ) = calpha*d( i )
799 50 CONTINUE
800*
801 END IF
802*
803* If matrix Hermitian, make D real
804*
805 IF( isym.EQ.0 ) THEN
806 DO 60 i = 1, mnmin
807 d( i ) = real( d( i ) )
808 60 CONTINUE
809 END IF
810*
811* Compute DL if grading set
812*
813 IF( igrade.EQ.1 .OR. igrade.EQ.3 .OR. igrade.EQ.4 .OR. igrade.EQ.
814 $ 5 .OR. igrade.EQ.6 ) THEN
815 CALL clatm1( model, condl, 0, idist, iseed, dl, m, info )
816 IF( info.NE.0 ) THEN
817 info = 3
818 RETURN
819 END IF
820 END IF
821*
822* Compute DR if grading set
823*
824 IF( igrade.EQ.2 .OR. igrade.EQ.3 ) THEN
825 CALL clatm1( moder, condr, 0, idist, iseed, dr, n, info )
826 IF( info.NE.0 ) THEN
827 info = 4
828 RETURN
829 END IF
830 END IF
831*
832* 3) Generate IWORK if pivoting
833*
834 IF( ipvtng.GT.0 ) THEN
835 DO 70 i = 1, npvts
836 iwork( i ) = i
837 70 CONTINUE
838 IF( fulbnd ) THEN
839 DO 80 i = 1, npvts
840 k = ipivot( i )
841 j = iwork( i )
842 iwork( i ) = iwork( k )
843 iwork( k ) = j
844 80 CONTINUE
845 ELSE
846 DO 90 i = npvts, 1, -1
847 k = ipivot( i )
848 j = iwork( i )
849 iwork( i ) = iwork( k )
850 iwork( k ) = j
851 90 CONTINUE
852 END IF
853 END IF
854*
855* 4) Generate matrices for each kind of PACKing
856* Always sweep matrix columnwise (if symmetric, upper
857* half only) so that matrix generated does not depend
858* on PACK
859*
860 IF( fulbnd ) THEN
861*
862* Use CLATM3 so matrices generated with differing PIVOTing only
863* differ only in the order of their rows and/or columns.
864*
865 IF( ipack.EQ.0 ) THEN
866 IF( isym.EQ.0 ) THEN
867 DO 110 j = 1, n
868 DO 100 i = 1, j
869 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
870 $ idist, iseed, d, igrade, dl, dr, ipvtng,
871 $ iwork, sparse )
872 a( isub, jsub ) = ctemp
873 a( jsub, isub ) = conjg( ctemp )
874 100 CONTINUE
875 110 CONTINUE
876 ELSE IF( isym.EQ.1 ) THEN
877 DO 130 j = 1, n
878 DO 120 i = 1, m
879 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
880 $ idist, iseed, d, igrade, dl, dr, ipvtng,
881 $ iwork, sparse )
882 a( isub, jsub ) = ctemp
883 120 CONTINUE
884 130 CONTINUE
885 ELSE IF( isym.EQ.2 ) THEN
886 DO 150 j = 1, n
887 DO 140 i = 1, j
888 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
889 $ idist, iseed, d, igrade, dl, dr, ipvtng,
890 $ iwork, sparse )
891 a( isub, jsub ) = ctemp
892 a( jsub, isub ) = ctemp
893 140 CONTINUE
894 150 CONTINUE
895 END IF
896*
897 ELSE IF( ipack.EQ.1 ) THEN
898*
899 DO 170 j = 1, n
900 DO 160 i = 1, j
901 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
902 $ idist,
903 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
904 $ sparse )
905 mnsub = min( isub, jsub )
906 mxsub = max( isub, jsub )
907 IF( mxsub.EQ.isub .AND. isym.EQ.0 ) THEN
908 a( mnsub, mxsub ) = conjg( ctemp )
909 ELSE
910 a( mnsub, mxsub ) = ctemp
911 END IF
912 IF( mnsub.NE.mxsub )
913 $ a( mxsub, mnsub ) = czero
914 160 CONTINUE
915 170 CONTINUE
916*
917 ELSE IF( ipack.EQ.2 ) THEN
918*
919 DO 190 j = 1, n
920 DO 180 i = 1, j
921 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
922 $ idist,
923 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
924 $ sparse )
925 mnsub = min( isub, jsub )
926 mxsub = max( isub, jsub )
927 IF( mxsub.EQ.jsub .AND. isym.EQ.0 ) THEN
928 a( mxsub, mnsub ) = conjg( ctemp )
929 ELSE
930 a( mxsub, mnsub ) = ctemp
931 END IF
932 IF( mnsub.NE.mxsub )
933 $ a( mnsub, mxsub ) = czero
934 180 CONTINUE
935 190 CONTINUE
936*
937 ELSE IF( ipack.EQ.3 ) THEN
938*
939 DO 210 j = 1, n
940 DO 200 i = 1, j
941 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
942 $ idist,
943 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
944 $ sparse )
945*
946* Compute K = location of (ISUB,JSUB) entry in packed
947* array
948*
949 mnsub = min( isub, jsub )
950 mxsub = max( isub, jsub )
951 k = mxsub*( mxsub-1 ) / 2 + mnsub
952*
953* Convert K to (IISUB,JJSUB) location
954*
955 jjsub = ( k-1 ) / lda + 1
956 iisub = k - lda*( jjsub-1 )
957*
958 IF( mxsub.EQ.isub .AND. isym.EQ.0 ) THEN
959 a( iisub, jjsub ) = conjg( ctemp )
960 ELSE
961 a( iisub, jjsub ) = ctemp
962 END IF
963 200 CONTINUE
964 210 CONTINUE
965*
966 ELSE IF( ipack.EQ.4 ) THEN
967*
968 DO 230 j = 1, n
969 DO 220 i = 1, j
970 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
971 $ idist,
972 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
973 $ sparse )
974*
975* Compute K = location of (I,J) entry in packed array
976*
977 mnsub = min( isub, jsub )
978 mxsub = max( isub, jsub )
979 IF( mnsub.EQ.1 ) THEN
980 k = mxsub
981 ELSE
982 k = n*( n+1 ) / 2 - ( n-mnsub+1 )*( n-mnsub+2 ) /
983 $ 2 + mxsub - mnsub + 1
984 END IF
985*
986* Convert K to (IISUB,JJSUB) location
987*
988 jjsub = ( k-1 ) / lda + 1
989 iisub = k - lda*( jjsub-1 )
990*
991 IF( mxsub.EQ.jsub .AND. isym.EQ.0 ) THEN
992 a( iisub, jjsub ) = conjg( ctemp )
993 ELSE
994 a( iisub, jjsub ) = ctemp
995 END IF
996 220 CONTINUE
997 230 CONTINUE
998*
999 ELSE IF( ipack.EQ.5 ) THEN
1000*
1001 DO 250 j = 1, n
1002 DO 240 i = j - kuu, j
1003 IF( i.LT.1 ) THEN
1004 a( j-i+1, i+n ) = czero
1005 ELSE
1006 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
1007 $ idist, iseed, d, igrade, dl, dr, ipvtng,
1008 $ iwork, sparse )
1009 mnsub = min( isub, jsub )
1010 mxsub = max( isub, jsub )
1011 IF( mxsub.EQ.jsub .AND. isym.EQ.0 ) THEN
1012 a( mxsub-mnsub+1, mnsub ) = conjg( ctemp )
1013 ELSE
1014 a( mxsub-mnsub+1, mnsub ) = ctemp
1015 END IF
1016 END IF
1017 240 CONTINUE
1018 250 CONTINUE
1019*
1020 ELSE IF( ipack.EQ.6 ) THEN
1021*
1022 DO 270 j = 1, n
1023 DO 260 i = j - kuu, j
1024 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
1025 $ idist,
1026 $ iseed, d, igrade, dl, dr, ipvtng, iwork,
1027 $ sparse )
1028 mnsub = min( isub, jsub )
1029 mxsub = max( isub, jsub )
1030 IF( mxsub.EQ.isub .AND. isym.EQ.0 ) THEN
1031 a( mnsub-mxsub+kuu+1, mxsub ) = conjg( ctemp )
1032 ELSE
1033 a( mnsub-mxsub+kuu+1, mxsub ) = ctemp
1034 END IF
1035 260 CONTINUE
1036 270 CONTINUE
1037*
1038 ELSE IF( ipack.EQ.7 ) THEN
1039*
1040 IF( isym.NE.1 ) THEN
1041 DO 290 j = 1, n
1042 DO 280 i = j - kuu, j
1043 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
1044 $ idist, iseed, d, igrade, dl, dr, ipvtng,
1045 $ iwork, sparse )
1046 mnsub = min( isub, jsub )
1047 mxsub = max( isub, jsub )
1048 IF( i.LT.1 )
1049 $ a( j-i+1+kuu, i+n ) = czero
1050 IF( mxsub.EQ.isub .AND. isym.EQ.0 ) THEN
1051 a( mnsub-mxsub+kuu+1, mxsub ) = conjg( ctemp )
1052 ELSE
1053 a( mnsub-mxsub+kuu+1, mxsub ) = ctemp
1054 END IF
1055 IF( i.GE.1 .AND. mnsub.NE.mxsub ) THEN
1056 IF( mnsub.EQ.isub .AND. isym.EQ.0 ) THEN
1057 a( mxsub-mnsub+1+kuu,
1058 $ mnsub ) = conjg( ctemp )
1059 ELSE
1060 a( mxsub-mnsub+1+kuu, mnsub ) = ctemp
1061 END IF
1062 END IF
1063 280 CONTINUE
1064 290 CONTINUE
1065 ELSE IF( isym.EQ.1 ) THEN
1066 DO 310 j = 1, n
1067 DO 300 i = j - kuu, j + kll
1068 ctemp = clatm3( m, n, i, j, isub, jsub, kl, ku,
1069 $ idist, iseed, d, igrade, dl, dr, ipvtng,
1070 $ iwork, sparse )
1071 a( isub-jsub+kuu+1, jsub ) = ctemp
1072 300 CONTINUE
1073 310 CONTINUE
1074 END IF
1075*
1076 END IF
1077*
1078 ELSE
1079*
1080* Use CLATM2
1081*
1082 IF( ipack.EQ.0 ) THEN
1083 IF( isym.EQ.0 ) THEN
1084 DO 330 j = 1, n
1085 DO 320 i = 1, j
1086 a( i, j ) = clatm2( m, n, i, j, kl, ku, idist,
1087 $ iseed, d, igrade, dl, dr, ipvtng,
1088 $ iwork, sparse )
1089 a( j, i ) = conjg( a( i, j ) )
1090 320 CONTINUE
1091 330 CONTINUE
1092 ELSE IF( isym.EQ.1 ) THEN
1093 DO 350 j = 1, n
1094 DO 340 i = 1, m
1095 a( i, j ) = clatm2( m, n, i, j, kl, ku, idist,
1096 $ iseed, d, igrade, dl, dr, ipvtng,
1097 $ iwork, sparse )
1098 340 CONTINUE
1099 350 CONTINUE
1100 ELSE IF( isym.EQ.2 ) THEN
1101 DO 370 j = 1, n
1102 DO 360 i = 1, j
1103 a( i, j ) = clatm2( m, n, i, j, kl, ku, idist,
1104 $ iseed, d, igrade, dl, dr, ipvtng,
1105 $ iwork, sparse )
1106 a( j, i ) = a( i, j )
1107 360 CONTINUE
1108 370 CONTINUE
1109 END IF
1110*
1111 ELSE IF( ipack.EQ.1 ) THEN
1112*
1113 DO 390 j = 1, n
1114 DO 380 i = 1, j
1115 a( i, j ) = clatm2( m, n, i, j, kl, ku, idist,
1116 $ iseed,
1117 $ d, igrade, dl, dr, ipvtng, iwork, sparse )
1118 IF( i.NE.j )
1119 $ a( j, i ) = czero
1120 380 CONTINUE
1121 390 CONTINUE
1122*
1123 ELSE IF( ipack.EQ.2 ) THEN
1124*
1125 DO 410 j = 1, n
1126 DO 400 i = 1, j
1127 IF( isym.EQ.0 ) THEN
1128 a( j, i ) = conjg( clatm2( m, n, i, j, kl, ku,
1129 $ idist, iseed, d, igrade, dl, dr,
1130 $ ipvtng, iwork, sparse ) )
1131 ELSE
1132 a( j, i ) = clatm2( m, n, i, j, kl, ku, idist,
1133 $ iseed, d, igrade, dl, dr, ipvtng,
1134 $ iwork, sparse )
1135 END IF
1136 IF( i.NE.j )
1137 $ a( i, j ) = czero
1138 400 CONTINUE
1139 410 CONTINUE
1140*
1141 ELSE IF( ipack.EQ.3 ) THEN
1142*
1143 isub = 0
1144 jsub = 1
1145 DO 430 j = 1, n
1146 DO 420 i = 1, j
1147 isub = isub + 1
1148 IF( isub.GT.lda ) THEN
1149 isub = 1
1150 jsub = jsub + 1
1151 END IF
1152 a( isub, jsub ) = clatm2( m, n, i, j, kl, ku,
1153 $ idist,
1154 $ iseed, d, igrade, dl, dr, ipvtng,
1155 $ iwork, sparse )
1156 420 CONTINUE
1157 430 CONTINUE
1158*
1159 ELSE IF( ipack.EQ.4 ) THEN
1160*
1161 IF( isym.EQ.0 .OR. isym.EQ.2 ) THEN
1162 DO 450 j = 1, n
1163 DO 440 i = 1, j
1164*
1165* Compute K = location of (I,J) entry in packed array
1166*
1167 IF( i.EQ.1 ) THEN
1168 k = j
1169 ELSE
1170 k = n*( n+1 ) / 2 - ( n-i+1 )*( n-i+2 ) / 2 +
1171 $ j - i + 1
1172 END IF
1173*
1174* Convert K to (ISUB,JSUB) location
1175*
1176 jsub = ( k-1 ) / lda + 1
1177 isub = k - lda*( jsub-1 )
1178*
1179 a( isub, jsub ) = clatm2( m, n, i, j, kl, ku,
1180 $ idist, iseed, d, igrade, dl, dr,
1181 $ ipvtng, iwork, sparse )
1182 IF( isym.EQ.0 )
1183 $ a( isub, jsub ) = conjg( a( isub, jsub ) )
1184 440 CONTINUE
1185 450 CONTINUE
1186 ELSE
1187 isub = 0
1188 jsub = 1
1189 DO 470 j = 1, n
1190 DO 460 i = j, m
1191 isub = isub + 1
1192 IF( isub.GT.lda ) THEN
1193 isub = 1
1194 jsub = jsub + 1
1195 END IF
1196 a( isub, jsub ) = clatm2( m, n, i, j, kl, ku,
1197 $ idist, iseed, d, igrade, dl, dr,
1198 $ ipvtng, iwork, sparse )
1199 460 CONTINUE
1200 470 CONTINUE
1201 END IF
1202*
1203 ELSE IF( ipack.EQ.5 ) THEN
1204*
1205 DO 490 j = 1, n
1206 DO 480 i = j - kuu, j
1207 IF( i.LT.1 ) THEN
1208 a( j-i+1, i+n ) = czero
1209 ELSE
1210 IF( isym.EQ.0 ) THEN
1211 a( j-i+1, i ) = conjg( clatm2( m, n, i, j,
1212 $ kl,
1213 $ ku, idist, iseed, d, igrade, dl,
1214 $ dr, ipvtng, iwork, sparse ) )
1215 ELSE
1216 a( j-i+1, i ) = clatm2( m, n, i, j, kl, ku,
1217 $ idist, iseed, d, igrade, dl, dr,
1218 $ ipvtng, iwork, sparse )
1219 END IF
1220 END IF
1221 480 CONTINUE
1222 490 CONTINUE
1223*
1224 ELSE IF( ipack.EQ.6 ) THEN
1225*
1226 DO 510 j = 1, n
1227 DO 500 i = j - kuu, j
1228 a( i-j+kuu+1, j ) = clatm2( m, n, i, j, kl, ku,
1229 $ idist,
1230 $ iseed, d, igrade, dl, dr, ipvtng,
1231 $ iwork, sparse )
1232 500 CONTINUE
1233 510 CONTINUE
1234*
1235 ELSE IF( ipack.EQ.7 ) THEN
1236*
1237 IF( isym.NE.1 ) THEN
1238 DO 530 j = 1, n
1239 DO 520 i = j - kuu, j
1240 a( i-j+kuu+1, j ) = clatm2( m, n, i, j, kl, ku,
1241 $ idist, iseed, d, igrade, dl,
1242 $ dr, ipvtng, iwork, sparse )
1243 IF( i.LT.1 )
1244 $ a( j-i+1+kuu, i+n ) = czero
1245 IF( i.GE.1 .AND. i.NE.j ) THEN
1246 IF( isym.EQ.0 ) THEN
1247 a( j-i+1+kuu, i ) = conjg( a( i-j+kuu+1,
1248 $ j ) )
1249 ELSE
1250 a( j-i+1+kuu, i ) = a( i-j+kuu+1, j )
1251 END IF
1252 END IF
1253 520 CONTINUE
1254 530 CONTINUE
1255 ELSE IF( isym.EQ.1 ) THEN
1256 DO 550 j = 1, n
1257 DO 540 i = j - kuu, j + kll
1258 a( i-j+kuu+1, j ) = clatm2( m, n, i, j, kl, ku,
1259 $ idist, iseed, d, igrade, dl,
1260 $ dr, ipvtng, iwork, sparse )
1261 540 CONTINUE
1262 550 CONTINUE
1263 END IF
1264*
1265 END IF
1266*
1267 END IF
1268*
1269* 5) Scaling the norm
1270*
1271 IF( ipack.EQ.0 ) THEN
1272 onorm = clange( 'M', m, n, a, lda, tempa )
1273 ELSE IF( ipack.EQ.1 ) THEN
1274 onorm = clansy( 'M', 'U', n, a, lda, tempa )
1275 ELSE IF( ipack.EQ.2 ) THEN
1276 onorm = clansy( 'M', 'L', n, a, lda, tempa )
1277 ELSE IF( ipack.EQ.3 ) THEN
1278 onorm = clansp( 'M', 'U', n, a, tempa )
1279 ELSE IF( ipack.EQ.4 ) THEN
1280 onorm = clansp( 'M', 'L', n, a, tempa )
1281 ELSE IF( ipack.EQ.5 ) THEN
1282 onorm = clansb( 'M', 'L', n, kll, a, lda, tempa )
1283 ELSE IF( ipack.EQ.6 ) THEN
1284 onorm = clansb( 'M', 'U', n, kuu, a, lda, tempa )
1285 ELSE IF( ipack.EQ.7 ) THEN
1286 onorm = clangb( 'M', n, kll, kuu, a, lda, tempa )
1287 END IF
1288*
1289 IF( anorm.GE.zero ) THEN
1290*
1291 IF( anorm.GT.zero .AND. onorm.EQ.zero ) THEN
1292*
1293* Desired scaling impossible
1294*
1295 info = 5
1296 RETURN
1297*
1298 ELSE IF( ( anorm.GT.one .AND. onorm.LT.one ) .OR.
1299 $ ( anorm.LT.one .AND. onorm.GT.one ) ) THEN
1300*
1301* Scale carefully to avoid over / underflow
1302*
1303 IF( ipack.LE.2 ) THEN
1304 DO 560 j = 1, n
1305 CALL csscal( m, one / onorm, a( 1, j ), 1 )
1306 CALL csscal( m, anorm, a( 1, j ), 1 )
1307 560 CONTINUE
1308*
1309 ELSE IF( ipack.EQ.3 .OR. ipack.EQ.4 ) THEN
1310*
1311 CALL csscal( n*( n+1 ) / 2, one / onorm, a, 1 )
1312 CALL csscal( n*( n+1 ) / 2, anorm, a, 1 )
1313*
1314 ELSE IF( ipack.GE.5 ) THEN
1315*
1316 DO 570 j = 1, n
1317 CALL csscal( kll+kuu+1, one / onorm, a( 1, j ), 1 )
1318 CALL csscal( kll+kuu+1, anorm, a( 1, j ), 1 )
1319 570 CONTINUE
1320*
1321 END IF
1322*
1323 ELSE
1324*
1325* Scale straightforwardly
1326*
1327 IF( ipack.LE.2 ) THEN
1328 DO 580 j = 1, n
1329 CALL csscal( m, anorm / onorm, a( 1, j ), 1 )
1330 580 CONTINUE
1331*
1332 ELSE IF( ipack.EQ.3 .OR. ipack.EQ.4 ) THEN
1333*
1334 CALL csscal( n*( n+1 ) / 2, anorm / onorm, a, 1 )
1335*
1336 ELSE IF( ipack.GE.5 ) THEN
1337*
1338 DO 590 j = 1, n
1339 CALL csscal( kll+kuu+1, anorm / onorm, a( 1, j ), 1 )
1340 590 CONTINUE
1341 END IF
1342*
1343 END IF
1344*
1345 END IF
1346*
1347* End of CLATMR
1348*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine clatm1(mode, cond, irsign, idist, iseed, d, n, info)
CLATM1
Definition clatm1.f:138
complex function clatm2(m, n, i, j, kl, ku, idist, iseed, d, igrade, dl, dr, ipvtng, iwork, sparse)
CLATM2
Definition clatm2.f:212
complex function clatm3(m, n, i, j, isub, jsub, kl, ku, idist, iseed, d, igrade, dl, dr, ipvtng, iwork, sparse)
CLATM3
Definition clatm3.f:229
real function clangb(norm, n, kl, ku, ab, ldab, work)
CLANGB returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition clangb.f:123
real function clange(norm, m, n, a, lda, work)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition clange.f:113
real function clansb(norm, uplo, n, k, ab, ldab, work)
CLANSB returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition clansb.f:128
real function clansy(norm, uplo, n, a, lda, work)
CLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition clansy.f:121
real function clansp(norm, uplo, n, ap, work)
CLANSP returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition clansp.f:113
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
subroutine csscal(n, sa, cx, incx)
CSSCAL
Definition csscal.f:78
Here is the call graph for this function:
Here is the caller graph for this function: