File:  [local] / rpl / lapack / lapack / dlansp.f
Revision 1.19: download - view: text, annotated - select for diffs - revision graph
Mon Aug 7 08:38:55 2023 UTC (8 months, 3 weeks ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_35, rpl-4_1_34, HEAD
Première mise à jour de lapack et blas.

    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.
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at
    6: *            http://www.netlib.org/lapack/explore-html/
    7: *
    8: *> \htmlonly
    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">
   15: *> [TXT]</a>
   16: *> \endhtmlonly
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       DOUBLE PRECISION FUNCTION DLANSP( NORM, UPLO, N, AP, WORK )
   22: *
   23: *       .. Scalar Arguments ..
   24: *       CHARACTER          NORM, UPLO
   25: *       INTEGER            N
   26: *       ..
   27: *       .. Array Arguments ..
   28: *       DOUBLE PRECISION   AP( * ), WORK( * )
   29: *       ..
   30: *
   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: *
  105: *> \author Univ. of Tennessee
  106: *> \author Univ. of California Berkeley
  107: *> \author Univ. of Colorado Denver
  108: *> \author NAG Ltd.
  109: *
  110: *> \ingroup doubleOTHERauxiliary
  111: *
  112: *  =====================================================================
  113:       DOUBLE PRECISION FUNCTION DLANSP( NORM, UPLO, N, AP, WORK )
  114: *
  115: *  -- LAPACK auxiliary routine --
  116: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  117: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  118: *
  119: *     .. Scalar Arguments ..
  120:       CHARACTER          NORM, UPLO
  121:       INTEGER            N
  122: *     ..
  123: *     .. Array Arguments ..
  124:       DOUBLE PRECISION   AP( * ), WORK( * )
  125: *     ..
  126: *
  127: * =====================================================================
  128: *
  129: *     .. Parameters ..
  130:       DOUBLE PRECISION   ONE, ZERO
  131:       PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
  132: *     ..
  133: *     .. Local Scalars ..
  134:       INTEGER            I, J, K
  135:       DOUBLE PRECISION   ABSA, SCALE, SUM, VALUE
  136: *     ..
  137: *     .. External Subroutines ..
  138:       EXTERNAL           DLASSQ
  139: *     ..
  140: *     .. External Functions ..
  141:       LOGICAL            LSAME, DISNAN
  142:       EXTERNAL           LSAME, DISNAN
  143: *     ..
  144: *     .. Intrinsic Functions ..
  145:       INTRINSIC          ABS, SQRT
  146: *     ..
  147: *     .. Executable Statements ..
  148: *
  149:       IF( N.EQ.0 ) THEN
  150:          VALUE = ZERO
  151:       ELSE IF( LSAME( NORM, 'M' ) ) THEN
  152: *
  153: *        Find max(abs(A(i,j))).
  154: *
  155:          VALUE = ZERO
  156:          IF( LSAME( UPLO, 'U' ) ) THEN
  157:             K = 1
  158:             DO 20 J = 1, N
  159:                DO 10 I = K, K + J - 1
  160:                   SUM = ABS( AP( I ) )
  161:                   IF( VALUE .LT. SUM .OR. DISNAN( SUM ) ) VALUE = SUM
  162:    10          CONTINUE
  163:                K = K + J
  164:    20       CONTINUE
  165:          ELSE
  166:             K = 1
  167:             DO 40 J = 1, N
  168:                DO 30 I = K, K + N - J
  169:                   SUM = ABS( AP( I ) )
  170:                   IF( VALUE .LT. SUM .OR. DISNAN( SUM ) ) VALUE = SUM
  171:    30          CONTINUE
  172:                K = K + N - J + 1
  173:    40       CONTINUE
  174:          END IF
  175:       ELSE IF( ( LSAME( NORM, 'I' ) ) .OR. ( LSAME( NORM, 'O' ) ) .OR.
  176:      $         ( NORM.EQ.'1' ) ) THEN
  177: *
  178: *        Find normI(A) ( = norm1(A), since A is symmetric).
  179: *
  180:          VALUE = ZERO
  181:          K = 1
  182:          IF( LSAME( UPLO, 'U' ) ) THEN
  183:             DO 60 J = 1, N
  184:                SUM = ZERO
  185:                DO 50 I = 1, J - 1
  186:                   ABSA = ABS( AP( K ) )
  187:                   SUM = SUM + ABSA
  188:                   WORK( I ) = WORK( I ) + ABSA
  189:                   K = K + 1
  190:    50          CONTINUE
  191:                WORK( J ) = SUM + ABS( AP( K ) )
  192:                K = K + 1
  193:    60       CONTINUE
  194:             DO 70 I = 1, N
  195:                SUM = WORK( I )
  196:                IF( VALUE .LT. SUM .OR. DISNAN( SUM ) ) VALUE = SUM
  197:    70       CONTINUE
  198:          ELSE
  199:             DO 80 I = 1, N
  200:                WORK( I ) = ZERO
  201:    80       CONTINUE
  202:             DO 100 J = 1, N
  203:                SUM = WORK( J ) + ABS( AP( K ) )
  204:                K = K + 1
  205:                DO 90 I = J + 1, N
  206:                   ABSA = ABS( AP( K ) )
  207:                   SUM = SUM + ABSA
  208:                   WORK( I ) = WORK( I ) + ABSA
  209:                   K = K + 1
  210:    90          CONTINUE
  211:                IF( VALUE .LT. SUM .OR. DISNAN( SUM ) ) VALUE = SUM
  212:   100       CONTINUE
  213:          END IF
  214:       ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
  215: *
  216: *        Find normF(A).
  217: *
  218:          SCALE = ZERO
  219:          SUM = ONE
  220:          K = 2
  221:          IF( LSAME( UPLO, 'U' ) ) THEN
  222:             DO 110 J = 2, N
  223:                CALL DLASSQ( J-1, AP( K ), 1, SCALE, SUM )
  224:                K = K + J
  225:   110       CONTINUE
  226:          ELSE
  227:             DO 120 J = 1, N - 1
  228:                CALL DLASSQ( N-J, AP( K ), 1, SCALE, SUM )
  229:                K = K + N - J + 1
  230:   120       CONTINUE
  231:          END IF
  232:          SUM = 2*SUM
  233:          K = 1
  234:          DO 130 I = 1, N
  235:             IF( AP( K ).NE.ZERO ) THEN
  236:                ABSA = ABS( AP( K ) )
  237:                IF( SCALE.LT.ABSA ) THEN
  238:                   SUM = ONE + SUM*( SCALE / ABSA )**2
  239:                   SCALE = ABSA
  240:                ELSE
  241:                   SUM = SUM + ( ABSA / SCALE )**2
  242:                END IF
  243:             END IF
  244:             IF( LSAME( UPLO, 'U' ) ) THEN
  245:                K = K + I + 1
  246:             ELSE
  247:                K = K + N - I + 1
  248:             END IF
  249:   130    CONTINUE
  250:          VALUE = SCALE*SQRT( SUM )
  251:       END IF
  252: *
  253:       DLANSP = VALUE
  254:       RETURN
  255: *
  256: *     End of DLANSP
  257: *
  258:       END

CVSweb interface <joel.bertrand@systella.fr>