Annotation of rpl/lapack/lapack/zlansp.f, revision 1.6

1.1       bertrand    1:       DOUBLE PRECISION FUNCTION ZLANSP( NORM, UPLO, N, AP, WORK )
                      2: *
                      3: *  -- LAPACK auxiliary routine (version 3.2) --
                      4: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                      5: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
                      6: *     November 2006
                      7: *
                      8: *     .. Scalar Arguments ..
                      9:       CHARACTER          NORM, UPLO
                     10:       INTEGER            N
                     11: *     ..
                     12: *     .. Array Arguments ..
                     13:       DOUBLE PRECISION   WORK( * )
                     14:       COMPLEX*16         AP( * )
                     15: *     ..
                     16: *
                     17: *  Purpose
                     18: *  =======
                     19: *
                     20: *  ZLANSP  returns the value of the one norm,  or the Frobenius norm, or
                     21: *  the  infinity norm,  or the  element of  largest absolute value  of a
                     22: *  complex symmetric matrix A,  supplied in packed form.
                     23: *
                     24: *  Description
                     25: *  ===========
                     26: *
                     27: *  ZLANSP returns the value
                     28: *
                     29: *     ZLANSP = ( max(abs(A(i,j))), NORM = 'M' or 'm'
                     30: *              (
                     31: *              ( norm1(A),         NORM = '1', 'O' or 'o'
                     32: *              (
                     33: *              ( normI(A),         NORM = 'I' or 'i'
                     34: *              (
                     35: *              ( normF(A),         NORM = 'F', 'f', 'E' or 'e'
                     36: *
                     37: *  where  norm1  denotes the  one norm of a matrix (maximum column sum),
                     38: *  normI  denotes the  infinity norm  of a matrix  (maximum row sum) and
                     39: *  normF  denotes the  Frobenius norm of a matrix (square root of sum of
                     40: *  squares).  Note that  max(abs(A(i,j)))  is not a consistent matrix norm.
                     41: *
                     42: *  Arguments
                     43: *  =========
                     44: *
                     45: *  NORM    (input) CHARACTER*1
                     46: *          Specifies the value to be returned in ZLANSP as described
                     47: *          above.
                     48: *
                     49: *  UPLO    (input) CHARACTER*1
                     50: *          Specifies whether the upper or lower triangular part of the
                     51: *          symmetric matrix A is supplied.
                     52: *          = 'U':  Upper triangular part of A is supplied
                     53: *          = 'L':  Lower triangular part of A is supplied
                     54: *
                     55: *  N       (input) INTEGER
                     56: *          The order of the matrix A.  N >= 0.  When N = 0, ZLANSP is
                     57: *          set to zero.
                     58: *
                     59: *  AP      (input) COMPLEX*16 array, dimension (N*(N+1)/2)
                     60: *          The upper or lower triangle of the symmetric matrix A, packed
                     61: *          columnwise in a linear array.  The j-th column of A is stored
                     62: *          in the array AP as follows:
                     63: *          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
                     64: *          if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n.
                     65: *
                     66: *  WORK    (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK)),
                     67: *          where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise,
                     68: *          WORK is not referenced.
                     69: *
                     70: * =====================================================================
                     71: *
                     72: *     .. Parameters ..
                     73:       DOUBLE PRECISION   ONE, ZERO
                     74:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
                     75: *     ..
                     76: *     .. Local Scalars ..
                     77:       INTEGER            I, J, K
                     78:       DOUBLE PRECISION   ABSA, SCALE, SUM, VALUE
                     79: *     ..
                     80: *     .. External Functions ..
                     81:       LOGICAL            LSAME
                     82:       EXTERNAL           LSAME
                     83: *     ..
                     84: *     .. External Subroutines ..
                     85:       EXTERNAL           ZLASSQ
                     86: *     ..
                     87: *     .. Intrinsic Functions ..
                     88:       INTRINSIC          ABS, DBLE, DIMAG, MAX, SQRT
                     89: *     ..
                     90: *     .. Executable Statements ..
                     91: *
                     92:       IF( N.EQ.0 ) THEN
                     93:          VALUE = ZERO
                     94:       ELSE IF( LSAME( NORM, 'M' ) ) THEN
                     95: *
                     96: *        Find max(abs(A(i,j))).
                     97: *
                     98:          VALUE = ZERO
                     99:          IF( LSAME( UPLO, 'U' ) ) THEN
                    100:             K = 1
                    101:             DO 20 J = 1, N
                    102:                DO 10 I = K, K + J - 1
                    103:                   VALUE = MAX( VALUE, ABS( AP( I ) ) )
                    104:    10          CONTINUE
                    105:                K = K + J
                    106:    20       CONTINUE
                    107:          ELSE
                    108:             K = 1
                    109:             DO 40 J = 1, N
                    110:                DO 30 I = K, K + N - J
                    111:                   VALUE = MAX( VALUE, ABS( AP( I ) ) )
                    112:    30          CONTINUE
                    113:                K = K + N - J + 1
                    114:    40       CONTINUE
                    115:          END IF
                    116:       ELSE IF( ( LSAME( NORM, 'I' ) ) .OR. ( LSAME( NORM, 'O' ) ) .OR.
                    117:      $         ( NORM.EQ.'1' ) ) THEN
                    118: *
                    119: *        Find normI(A) ( = norm1(A), since A is symmetric).
                    120: *
                    121:          VALUE = ZERO
                    122:          K = 1
                    123:          IF( LSAME( UPLO, 'U' ) ) THEN
                    124:             DO 60 J = 1, N
                    125:                SUM = ZERO
                    126:                DO 50 I = 1, J - 1
                    127:                   ABSA = ABS( AP( K ) )
                    128:                   SUM = SUM + ABSA
                    129:                   WORK( I ) = WORK( I ) + ABSA
                    130:                   K = K + 1
                    131:    50          CONTINUE
                    132:                WORK( J ) = SUM + ABS( AP( K ) )
                    133:                K = K + 1
                    134:    60       CONTINUE
                    135:             DO 70 I = 1, N
                    136:                VALUE = MAX( VALUE, WORK( I ) )
                    137:    70       CONTINUE
                    138:          ELSE
                    139:             DO 80 I = 1, N
                    140:                WORK( I ) = ZERO
                    141:    80       CONTINUE
                    142:             DO 100 J = 1, N
                    143:                SUM = WORK( J ) + ABS( AP( K ) )
                    144:                K = K + 1
                    145:                DO 90 I = J + 1, N
                    146:                   ABSA = ABS( AP( K ) )
                    147:                   SUM = SUM + ABSA
                    148:                   WORK( I ) = WORK( I ) + ABSA
                    149:                   K = K + 1
                    150:    90          CONTINUE
                    151:                VALUE = MAX( VALUE, SUM )
                    152:   100       CONTINUE
                    153:          END IF
                    154:       ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
                    155: *
                    156: *        Find normF(A).
                    157: *
                    158:          SCALE = ZERO
                    159:          SUM = ONE
                    160:          K = 2
                    161:          IF( LSAME( UPLO, 'U' ) ) THEN
                    162:             DO 110 J = 2, N
                    163:                CALL ZLASSQ( J-1, AP( K ), 1, SCALE, SUM )
                    164:                K = K + J
                    165:   110       CONTINUE
                    166:          ELSE
                    167:             DO 120 J = 1, N - 1
                    168:                CALL ZLASSQ( N-J, AP( K ), 1, SCALE, SUM )
                    169:                K = K + N - J + 1
                    170:   120       CONTINUE
                    171:          END IF
                    172:          SUM = 2*SUM
                    173:          K = 1
                    174:          DO 130 I = 1, N
                    175:             IF( DBLE( AP( K ) ).NE.ZERO ) THEN
                    176:                ABSA = ABS( DBLE( AP( K ) ) )
                    177:                IF( SCALE.LT.ABSA ) THEN
                    178:                   SUM = ONE + SUM*( SCALE / ABSA )**2
                    179:                   SCALE = ABSA
                    180:                ELSE
                    181:                   SUM = SUM + ( ABSA / SCALE )**2
                    182:                END IF
                    183:             END IF
                    184:             IF( DIMAG( AP( K ) ).NE.ZERO ) THEN
                    185:                ABSA = ABS( DIMAG( AP( K ) ) )
                    186:                IF( SCALE.LT.ABSA ) THEN
                    187:                   SUM = ONE + SUM*( SCALE / ABSA )**2
                    188:                   SCALE = ABSA
                    189:                ELSE
                    190:                   SUM = SUM + ( ABSA / SCALE )**2
                    191:                END IF
                    192:             END IF
                    193:             IF( LSAME( UPLO, 'U' ) ) THEN
                    194:                K = K + I + 1
                    195:             ELSE
                    196:                K = K + N - I + 1
                    197:             END IF
                    198:   130    CONTINUE
                    199:          VALUE = SCALE*SQRT( SUM )
                    200:       END IF
                    201: *
                    202:       ZLANSP = VALUE
                    203:       RETURN
                    204: *
                    205: *     End of ZLANSP
                    206: *
                    207:       END

CVSweb interface <joel.bertrand@systella.fr>