Annotation of rpl/lapack/lapack/dlacon.f, revision 1.19

1.12      bertrand    1: *> \brief \b DLACON estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vector products.
1.9       bertrand    2: *
                      3: *  =========== DOCUMENTATION ===========
                      4: *
1.16      bertrand    5: * Online html documentation available at
                      6: *            http://www.netlib.org/lapack/explore-html/
1.9       bertrand    7: *
                      8: *> \htmlonly
1.16      bertrand    9: *> Download DLACON + dependencies
                     10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlacon.f">
                     11: *> [TGZ]</a>
                     12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlacon.f">
                     13: *> [ZIP]</a>
                     14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlacon.f">
1.9       bertrand   15: *> [TXT]</a>
1.16      bertrand   16: *> \endhtmlonly
1.9       bertrand   17: *
                     18: *  Definition:
                     19: *  ===========
                     20: *
                     21: *       SUBROUTINE DLACON( N, V, X, ISGN, EST, KASE )
1.16      bertrand   22: *
1.9       bertrand   23: *       .. Scalar Arguments ..
                     24: *       INTEGER            KASE, N
                     25: *       DOUBLE PRECISION   EST
                     26: *       ..
                     27: *       .. Array Arguments ..
                     28: *       INTEGER            ISGN( * )
                     29: *       DOUBLE PRECISION   V( * ), X( * )
                     30: *       ..
1.16      bertrand   31: *
1.9       bertrand   32: *
                     33: *> \par Purpose:
                     34: *  =============
                     35: *>
                     36: *> \verbatim
                     37: *>
                     38: *> DLACON estimates the 1-norm of a square, real matrix A.
                     39: *> Reverse communication is used for evaluating matrix-vector products.
                     40: *> \endverbatim
                     41: *
                     42: *  Arguments:
                     43: *  ==========
                     44: *
                     45: *> \param[in] N
                     46: *> \verbatim
                     47: *>          N is INTEGER
                     48: *>         The order of the matrix.  N >= 1.
                     49: *> \endverbatim
                     50: *>
                     51: *> \param[out] V
                     52: *> \verbatim
                     53: *>          V is DOUBLE PRECISION array, dimension (N)
                     54: *>         On the final return, V = A*W,  where  EST = norm(V)/norm(W)
                     55: *>         (W is not returned).
                     56: *> \endverbatim
                     57: *>
                     58: *> \param[in,out] X
                     59: *> \verbatim
                     60: *>          X is DOUBLE PRECISION array, dimension (N)
                     61: *>         On an intermediate return, X should be overwritten by
                     62: *>               A * X,   if KASE=1,
                     63: *>               A**T * X,  if KASE=2,
                     64: *>         and DLACON must be re-called with all the other parameters
                     65: *>         unchanged.
                     66: *> \endverbatim
                     67: *>
                     68: *> \param[out] ISGN
                     69: *> \verbatim
                     70: *>          ISGN is INTEGER array, dimension (N)
                     71: *> \endverbatim
                     72: *>
                     73: *> \param[in,out] EST
                     74: *> \verbatim
                     75: *>          EST is DOUBLE PRECISION
                     76: *>         On entry with KASE = 1 or 2 and JUMP = 3, EST should be
                     77: *>         unchanged from the previous call to DLACON.
1.16      bertrand   78: *>         On exit, EST is an estimate (a lower bound) for norm(A).
1.9       bertrand   79: *> \endverbatim
                     80: *>
                     81: *> \param[in,out] KASE
                     82: *> \verbatim
                     83: *>          KASE is INTEGER
                     84: *>         On the initial call to DLACON, KASE should be 0.
                     85: *>         On an intermediate return, KASE will be 1 or 2, indicating
                     86: *>         whether X should be overwritten by A * X  or A**T * X.
                     87: *>         On the final return from DLACON, KASE will again be 0.
                     88: *> \endverbatim
                     89: *
                     90: *  Authors:
                     91: *  ========
                     92: *
1.16      bertrand   93: *> \author Univ. of Tennessee
                     94: *> \author Univ. of California Berkeley
                     95: *> \author Univ. of Colorado Denver
                     96: *> \author NAG Ltd.
1.9       bertrand   97: *
                     98: *> \ingroup doubleOTHERauxiliary
                     99: *
                    100: *> \par Contributors:
                    101: *  ==================
                    102: *>
                    103: *>  Nick Higham, University of Manchester. \n
                    104: *>  Originally named SONEST, dated March 16, 1988.
                    105: *
                    106: *> \par References:
                    107: *  ================
                    108: *>
                    109: *>  N.J. Higham, "FORTRAN codes for estimating the one-norm of
                    110: *>  a real or complex matrix, with applications to condition estimation",
                    111: *>  ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988.
                    112: *>
                    113: *  =====================================================================
1.1       bertrand  114:       SUBROUTINE DLACON( N, V, X, ISGN, EST, KASE )
                    115: *
1.19    ! bertrand  116: *  -- LAPACK auxiliary routine --
1.1       bertrand  117: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                    118: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
                    119: *
                    120: *     .. Scalar Arguments ..
                    121:       INTEGER            KASE, N
                    122:       DOUBLE PRECISION   EST
                    123: *     ..
                    124: *     .. Array Arguments ..
                    125:       INTEGER            ISGN( * )
                    126:       DOUBLE PRECISION   V( * ), X( * )
                    127: *     ..
                    128: *
                    129: *  =====================================================================
                    130: *
                    131: *     .. Parameters ..
                    132:       INTEGER            ITMAX
                    133:       PARAMETER          ( ITMAX = 5 )
                    134:       DOUBLE PRECISION   ZERO, ONE, TWO
                    135:       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0, TWO = 2.0D+0 )
                    136: *     ..
                    137: *     .. Local Scalars ..
                    138:       INTEGER            I, ITER, J, JLAST, JUMP
                    139:       DOUBLE PRECISION   ALTSGN, ESTOLD, TEMP
                    140: *     ..
                    141: *     .. External Functions ..
                    142:       INTEGER            IDAMAX
                    143:       DOUBLE PRECISION   DASUM
                    144:       EXTERNAL           IDAMAX, DASUM
                    145: *     ..
                    146: *     .. External Subroutines ..
                    147:       EXTERNAL           DCOPY
                    148: *     ..
                    149: *     .. Intrinsic Functions ..
                    150:       INTRINSIC          ABS, DBLE, NINT, SIGN
                    151: *     ..
                    152: *     .. Save statement ..
                    153:       SAVE
                    154: *     ..
                    155: *     .. Executable Statements ..
                    156: *
                    157:       IF( KASE.EQ.0 ) THEN
                    158:          DO 10 I = 1, N
                    159:             X( I ) = ONE / DBLE( N )
                    160:    10    CONTINUE
                    161:          KASE = 1
                    162:          JUMP = 1
                    163:          RETURN
                    164:       END IF
                    165: *
                    166:       GO TO ( 20, 40, 70, 110, 140 )JUMP
                    167: *
                    168: *     ................ ENTRY   (JUMP = 1)
                    169: *     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY A*X.
                    170: *
                    171:    20 CONTINUE
                    172:       IF( N.EQ.1 ) THEN
                    173:          V( 1 ) = X( 1 )
                    174:          EST = ABS( V( 1 ) )
                    175: *        ... QUIT
                    176:          GO TO 150
                    177:       END IF
                    178:       EST = DASUM( N, X, 1 )
                    179: *
                    180:       DO 30 I = 1, N
                    181:          X( I ) = SIGN( ONE, X( I ) )
                    182:          ISGN( I ) = NINT( X( I ) )
                    183:    30 CONTINUE
                    184:       KASE = 2
                    185:       JUMP = 2
                    186:       RETURN
                    187: *
                    188: *     ................ ENTRY   (JUMP = 2)
                    189: *     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY TRANSPOSE(A)*X.
                    190: *
                    191:    40 CONTINUE
                    192:       J = IDAMAX( N, X, 1 )
                    193:       ITER = 2
                    194: *
                    195: *     MAIN LOOP - ITERATIONS 2,3,...,ITMAX.
                    196: *
                    197:    50 CONTINUE
                    198:       DO 60 I = 1, N
                    199:          X( I ) = ZERO
                    200:    60 CONTINUE
                    201:       X( J ) = ONE
                    202:       KASE = 1
                    203:       JUMP = 3
                    204:       RETURN
                    205: *
                    206: *     ................ ENTRY   (JUMP = 3)
                    207: *     X HAS BEEN OVERWRITTEN BY A*X.
                    208: *
                    209:    70 CONTINUE
                    210:       CALL DCOPY( N, X, 1, V, 1 )
                    211:       ESTOLD = EST
                    212:       EST = DASUM( N, V, 1 )
                    213:       DO 80 I = 1, N
                    214:          IF( NINT( SIGN( ONE, X( I ) ) ).NE.ISGN( I ) )
                    215:      $      GO TO 90
                    216:    80 CONTINUE
                    217: *     REPEATED SIGN VECTOR DETECTED, HENCE ALGORITHM HAS CONVERGED.
                    218:       GO TO 120
                    219: *
                    220:    90 CONTINUE
                    221: *     TEST FOR CYCLING.
                    222:       IF( EST.LE.ESTOLD )
                    223:      $   GO TO 120
                    224: *
                    225:       DO 100 I = 1, N
                    226:          X( I ) = SIGN( ONE, X( I ) )
                    227:          ISGN( I ) = NINT( X( I ) )
                    228:   100 CONTINUE
                    229:       KASE = 2
                    230:       JUMP = 4
                    231:       RETURN
                    232: *
                    233: *     ................ ENTRY   (JUMP = 4)
                    234: *     X HAS BEEN OVERWRITTEN BY TRANSPOSE(A)*X.
                    235: *
                    236:   110 CONTINUE
                    237:       JLAST = J
                    238:       J = IDAMAX( N, X, 1 )
                    239:       IF( ( X( JLAST ).NE.ABS( X( J ) ) ) .AND. ( ITER.LT.ITMAX ) ) THEN
                    240:          ITER = ITER + 1
                    241:          GO TO 50
                    242:       END IF
                    243: *
                    244: *     ITERATION COMPLETE.  FINAL STAGE.
                    245: *
                    246:   120 CONTINUE
                    247:       ALTSGN = ONE
                    248:       DO 130 I = 1, N
                    249:          X( I ) = ALTSGN*( ONE+DBLE( I-1 ) / DBLE( N-1 ) )
                    250:          ALTSGN = -ALTSGN
                    251:   130 CONTINUE
                    252:       KASE = 1
                    253:       JUMP = 5
                    254:       RETURN
                    255: *
                    256: *     ................ ENTRY   (JUMP = 5)
                    257: *     X HAS BEEN OVERWRITTEN BY A*X.
                    258: *
                    259:   140 CONTINUE
                    260:       TEMP = TWO*( DASUM( N, X, 1 ) / DBLE( 3*N ) )
                    261:       IF( TEMP.GT.EST ) THEN
                    262:          CALL DCOPY( N, X, 1, V, 1 )
                    263:          EST = TEMP
                    264:       END IF
                    265: *
                    266:   150 CONTINUE
                    267:       KASE = 0
                    268:       RETURN
                    269: *
                    270: *     End of DLACON
                    271: *
                    272:       END

CVSweb interface <joel.bertrand@systella.fr>