Annotation of rpl/lapack/lapack/dlansp.f, revision 1.18

1.11      bertrand    1: *> \brief \b DLANSP returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a symmetric matrix supplied in packed form.
1.8       bertrand    2: *
                      3: *  =========== DOCUMENTATION ===========
                      4: *
1.15      bertrand    5: * Online html documentation available at
                      6: *            http://www.netlib.org/lapack/explore-html/
1.8       bertrand    7: *
                      8: *> \htmlonly
1.15      bertrand    9: *> Download DLANSP + dependencies
                     10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlansp.f">
                     11: *> [TGZ]</a>
                     12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlansp.f">
                     13: *> [ZIP]</a>
                     14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlansp.f">
1.8       bertrand   15: *> [TXT]</a>
1.15      bertrand   16: *> \endhtmlonly
1.8       bertrand   17: *
                     18: *  Definition:
                     19: *  ===========
                     20: *
                     21: *       DOUBLE PRECISION FUNCTION DLANSP( NORM, UPLO, N, AP, WORK )
1.15      bertrand   22: *
1.8       bertrand   23: *       .. Scalar Arguments ..
                     24: *       CHARACTER          NORM, UPLO
                     25: *       INTEGER            N
                     26: *       ..
                     27: *       .. Array Arguments ..
                     28: *       DOUBLE PRECISION   AP( * ), WORK( * )
                     29: *       ..
1.15      bertrand   30: *
1.8       bertrand   31: *
                     32: *> \par Purpose:
                     33: *  =============
                     34: *>
                     35: *> \verbatim
                     36: *>
                     37: *> DLANSP  returns the value of the one norm,  or the Frobenius norm, or
                     38: *> the  infinity norm,  or the  element of  largest absolute value  of a
                     39: *> real symmetric matrix A,  supplied in packed form.
                     40: *> \endverbatim
                     41: *>
                     42: *> \return DLANSP
                     43: *> \verbatim
                     44: *>
                     45: *>    DLANSP = ( max(abs(A(i,j))), NORM = 'M' or 'm'
                     46: *>             (
                     47: *>             ( norm1(A),         NORM = '1', 'O' or 'o'
                     48: *>             (
                     49: *>             ( normI(A),         NORM = 'I' or 'i'
                     50: *>             (
                     51: *>             ( normF(A),         NORM = 'F', 'f', 'E' or 'e'
                     52: *>
                     53: *> where  norm1  denotes the  one norm of a matrix (maximum column sum),
                     54: *> normI  denotes the  infinity norm  of a matrix  (maximum row sum) and
                     55: *> normF  denotes the  Frobenius norm of a matrix (square root of sum of
                     56: *> squares).  Note that  max(abs(A(i,j)))  is not a consistent matrix norm.
                     57: *> \endverbatim
                     58: *
                     59: *  Arguments:
                     60: *  ==========
                     61: *
                     62: *> \param[in] NORM
                     63: *> \verbatim
                     64: *>          NORM is CHARACTER*1
                     65: *>          Specifies the value to be returned in DLANSP as described
                     66: *>          above.
                     67: *> \endverbatim
                     68: *>
                     69: *> \param[in] UPLO
                     70: *> \verbatim
                     71: *>          UPLO is CHARACTER*1
                     72: *>          Specifies whether the upper or lower triangular part of the
                     73: *>          symmetric matrix A is supplied.
                     74: *>          = 'U':  Upper triangular part of A is supplied
                     75: *>          = 'L':  Lower triangular part of A is supplied
                     76: *> \endverbatim
                     77: *>
                     78: *> \param[in] N
                     79: *> \verbatim
                     80: *>          N is INTEGER
                     81: *>          The order of the matrix A.  N >= 0.  When N = 0, DLANSP is
                     82: *>          set to zero.
                     83: *> \endverbatim
                     84: *>
                     85: *> \param[in] AP
                     86: *> \verbatim
                     87: *>          AP is DOUBLE PRECISION array, dimension (N*(N+1)/2)
                     88: *>          The upper or lower triangle of the symmetric matrix A, packed
                     89: *>          columnwise in a linear array.  The j-th column of A is stored
                     90: *>          in the array AP as follows:
                     91: *>          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
                     92: *>          if UPLO = 'L', AP(i + (j-1)*(2n-j)/2) = A(i,j) for j<=i<=n.
                     93: *> \endverbatim
                     94: *>
                     95: *> \param[out] WORK
                     96: *> \verbatim
                     97: *>          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK)),
                     98: *>          where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise,
                     99: *>          WORK is not referenced.
                    100: *> \endverbatim
                    101: *
                    102: *  Authors:
                    103: *  ========
                    104: *
1.15      bertrand  105: *> \author Univ. of Tennessee
                    106: *> \author Univ. of California Berkeley
                    107: *> \author Univ. of Colorado Denver
                    108: *> \author NAG Ltd.
1.8       bertrand  109: *
1.15      bertrand  110: *> \date December 2016
1.8       bertrand  111: *
                    112: *> \ingroup doubleOTHERauxiliary
                    113: *
                    114: *  =====================================================================
1.1       bertrand  115:       DOUBLE PRECISION FUNCTION DLANSP( NORM, UPLO, N, AP, WORK )
                    116: *
1.15      bertrand  117: *  -- LAPACK auxiliary routine (version 3.7.0) --
1.1       bertrand  118: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                    119: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.15      bertrand  120: *     December 2016
1.1       bertrand  121: *
1.18    ! bertrand  122:       IMPLICIT NONE
1.1       bertrand  123: *     .. Scalar Arguments ..
                    124:       CHARACTER          NORM, UPLO
                    125:       INTEGER            N
                    126: *     ..
                    127: *     .. Array Arguments ..
                    128:       DOUBLE PRECISION   AP( * ), WORK( * )
                    129: *     ..
                    130: *
                    131: * =====================================================================
                    132: *
                    133: *     .. Parameters ..
                    134:       DOUBLE PRECISION   ONE, ZERO
                    135:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
                    136: *     ..
                    137: *     .. Local Scalars ..
                    138:       INTEGER            I, J, K
1.18    ! bertrand  139:       DOUBLE PRECISION   ABSA, SUM, VALUE
1.1       bertrand  140: *     ..
1.18    ! bertrand  141: *     .. Local Arrays ..
        !           142:       DOUBLE PRECISION   SSQ( 2 ), COLSSQ( 2 )
1.1       bertrand  143: *     ..
                    144: *     .. External Functions ..
1.11      bertrand  145:       LOGICAL            LSAME, DISNAN
                    146:       EXTERNAL           LSAME, DISNAN
1.1       bertrand  147: *     ..
1.18    ! bertrand  148: *     .. External Subroutines ..
        !           149:       EXTERNAL           DLASSQ, DCOMBSSQ
        !           150: *     ..
1.1       bertrand  151: *     .. Intrinsic Functions ..
1.11      bertrand  152:       INTRINSIC          ABS, SQRT
1.1       bertrand  153: *     ..
                    154: *     .. Executable Statements ..
                    155: *
                    156:       IF( N.EQ.0 ) THEN
                    157:          VALUE = ZERO
                    158:       ELSE IF( LSAME( NORM, 'M' ) ) THEN
                    159: *
                    160: *        Find max(abs(A(i,j))).
                    161: *
                    162:          VALUE = ZERO
                    163:          IF( LSAME( UPLO, 'U' ) ) THEN
                    164:             K = 1
                    165:             DO 20 J = 1, N
                    166:                DO 10 I = K, K + J - 1
1.11      bertrand  167:                   SUM = ABS( AP( I ) )
                    168:                   IF( VALUE .LT. SUM .OR. DISNAN( SUM ) ) VALUE = SUM
1.1       bertrand  169:    10          CONTINUE
                    170:                K = K + J
                    171:    20       CONTINUE
                    172:          ELSE
                    173:             K = 1
                    174:             DO 40 J = 1, N
                    175:                DO 30 I = K, K + N - J
1.11      bertrand  176:                   SUM = ABS( AP( I ) )
                    177:                   IF( VALUE .LT. SUM .OR. DISNAN( SUM ) ) VALUE = SUM
1.1       bertrand  178:    30          CONTINUE
                    179:                K = K + N - J + 1
                    180:    40       CONTINUE
                    181:          END IF
                    182:       ELSE IF( ( LSAME( NORM, 'I' ) ) .OR. ( LSAME( NORM, 'O' ) ) .OR.
                    183:      $         ( NORM.EQ.'1' ) ) THEN
                    184: *
                    185: *        Find normI(A) ( = norm1(A), since A is symmetric).
                    186: *
                    187:          VALUE = ZERO
                    188:          K = 1
                    189:          IF( LSAME( UPLO, 'U' ) ) THEN
                    190:             DO 60 J = 1, N
                    191:                SUM = ZERO
                    192:                DO 50 I = 1, J - 1
                    193:                   ABSA = ABS( AP( K ) )
                    194:                   SUM = SUM + ABSA
                    195:                   WORK( I ) = WORK( I ) + ABSA
                    196:                   K = K + 1
                    197:    50          CONTINUE
                    198:                WORK( J ) = SUM + ABS( AP( K ) )
                    199:                K = K + 1
                    200:    60       CONTINUE
                    201:             DO 70 I = 1, N
1.11      bertrand  202:                SUM = WORK( I )
                    203:                IF( VALUE .LT. SUM .OR. DISNAN( SUM ) ) VALUE = SUM
1.1       bertrand  204:    70       CONTINUE
                    205:          ELSE
                    206:             DO 80 I = 1, N
                    207:                WORK( I ) = ZERO
                    208:    80       CONTINUE
                    209:             DO 100 J = 1, N
                    210:                SUM = WORK( J ) + ABS( AP( K ) )
                    211:                K = K + 1
                    212:                DO 90 I = J + 1, N
                    213:                   ABSA = ABS( AP( K ) )
                    214:                   SUM = SUM + ABSA
                    215:                   WORK( I ) = WORK( I ) + ABSA
                    216:                   K = K + 1
                    217:    90          CONTINUE
1.11      bertrand  218:                IF( VALUE .LT. SUM .OR. DISNAN( SUM ) ) VALUE = SUM
1.1       bertrand  219:   100       CONTINUE
                    220:          END IF
                    221:       ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
                    222: *
                    223: *        Find normF(A).
1.18    ! bertrand  224: *        SSQ(1) is scale
        !           225: *        SSQ(2) is sum-of-squares
        !           226: *        For better accuracy, sum each column separately.
        !           227: *
        !           228:          SSQ( 1 ) = ZERO
        !           229:          SSQ( 2 ) = ONE
        !           230: *
        !           231: *        Sum off-diagonals
1.1       bertrand  232: *
                    233:          K = 2
                    234:          IF( LSAME( UPLO, 'U' ) ) THEN
                    235:             DO 110 J = 2, N
1.18    ! bertrand  236:                COLSSQ( 1 ) = ZERO
        !           237:                COLSSQ( 2 ) = ONE
        !           238:                CALL DLASSQ( J-1, AP( K ), 1, COLSSQ( 1 ), COLSSQ( 2 ) )
        !           239:                CALL DCOMBSSQ( SSQ, COLSSQ )
1.1       bertrand  240:                K = K + J
                    241:   110       CONTINUE
                    242:          ELSE
                    243:             DO 120 J = 1, N - 1
1.18    ! bertrand  244:                COLSSQ( 1 ) = ZERO
        !           245:                COLSSQ( 2 ) = ONE
        !           246:                CALL DLASSQ( N-J, AP( K ), 1, COLSSQ( 1 ), COLSSQ( 2 ) )
        !           247:                CALL DCOMBSSQ( SSQ, COLSSQ )
1.1       bertrand  248:                K = K + N - J + 1
                    249:   120       CONTINUE
                    250:          END IF
1.18    ! bertrand  251:          SSQ( 2 ) = 2*SSQ( 2 )
        !           252: *
        !           253: *        Sum diagonal
        !           254: *
1.1       bertrand  255:          K = 1
1.18    ! bertrand  256:          COLSSQ( 1 ) = ZERO
        !           257:          COLSSQ( 2 ) = ONE
1.1       bertrand  258:          DO 130 I = 1, N
                    259:             IF( AP( K ).NE.ZERO ) THEN
                    260:                ABSA = ABS( AP( K ) )
1.18    ! bertrand  261:                IF( COLSSQ( 1 ).LT.ABSA ) THEN
        !           262:                   COLSSQ( 2 ) = ONE + COLSSQ(2)*( COLSSQ(1) / ABSA )**2
        !           263:                   COLSSQ( 1 ) = ABSA
1.1       bertrand  264:                ELSE
1.18    ! bertrand  265:                   COLSSQ( 2 ) = COLSSQ( 2 ) + ( ABSA / COLSSQ( 1 ) )**2
1.1       bertrand  266:                END IF
                    267:             END IF
                    268:             IF( LSAME( UPLO, 'U' ) ) THEN
                    269:                K = K + I + 1
                    270:             ELSE
                    271:                K = K + N - I + 1
                    272:             END IF
                    273:   130    CONTINUE
1.18    ! bertrand  274:          CALL DCOMBSSQ( SSQ, COLSSQ )
        !           275:          VALUE = SSQ( 1 )*SQRT( SSQ( 2 ) )
1.1       bertrand  276:       END IF
                    277: *
                    278:       DLANSP = VALUE
                    279:       RETURN
                    280: *
                    281: *     End of DLANSP
                    282: *
                    283:       END

CVSweb interface <joel.bertrand@systella.fr>