File:  [local] / rpl / lapack / lapack / dlacon.f
Revision 1.15: download - view: text, annotated - select for diffs - revision graph
Sat Aug 27 15:34:25 2016 UTC (7 years, 9 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_25, HEAD
Cohérence Lapack.

    1: *> \brief \b DLACON estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vector products.
    2: *
    3: *  =========== DOCUMENTATION ===========
    4: *
    5: * Online html documentation available at 
    6: *            http://www.netlib.org/lapack/explore-html/ 
    7: *
    8: *> \htmlonly
    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"> 
   15: *> [TXT]</a>
   16: *> \endhtmlonly 
   17: *
   18: *  Definition:
   19: *  ===========
   20: *
   21: *       SUBROUTINE DLACON( N, V, X, ISGN, EST, KASE )
   22:    23: *       .. Scalar Arguments ..
   24: *       INTEGER            KASE, N
   25: *       DOUBLE PRECISION   EST
   26: *       ..
   27: *       .. Array Arguments ..
   28: *       INTEGER            ISGN( * )
   29: *       DOUBLE PRECISION   V( * ), X( * )
   30: *       ..
   31: *  
   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.
   78: *>         On exit, EST is an estimate (a lower bound) for norm(A). 
   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: *
   93: *> \author Univ. of Tennessee 
   94: *> \author Univ. of California Berkeley 
   95: *> \author Univ. of Colorado Denver 
   96: *> \author NAG Ltd. 
   97: *
   98: *> \date September 2012
   99: *
  100: *> \ingroup doubleOTHERauxiliary
  101: *
  102: *> \par Contributors:
  103: *  ==================
  104: *>
  105: *>  Nick Higham, University of Manchester. \n
  106: *>  Originally named SONEST, dated March 16, 1988.
  107: *
  108: *> \par References:
  109: *  ================
  110: *>
  111: *>  N.J. Higham, "FORTRAN codes for estimating the one-norm of
  112: *>  a real or complex matrix, with applications to condition estimation",
  113: *>  ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988.
  114: *>
  115: *  =====================================================================
  116:       SUBROUTINE DLACON( N, V, X, ISGN, EST, KASE )
  117: *
  118: *  -- LAPACK auxiliary routine (version 3.4.2) --
  119: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  120: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  121: *     September 2012
  122: *
  123: *     .. Scalar Arguments ..
  124:       INTEGER            KASE, N
  125:       DOUBLE PRECISION   EST
  126: *     ..
  127: *     .. Array Arguments ..
  128:       INTEGER            ISGN( * )
  129:       DOUBLE PRECISION   V( * ), X( * )
  130: *     ..
  131: *
  132: *  =====================================================================
  133: *
  134: *     .. Parameters ..
  135:       INTEGER            ITMAX
  136:       PARAMETER          ( ITMAX = 5 )
  137:       DOUBLE PRECISION   ZERO, ONE, TWO
  138:       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0, TWO = 2.0D+0 )
  139: *     ..
  140: *     .. Local Scalars ..
  141:       INTEGER            I, ITER, J, JLAST, JUMP
  142:       DOUBLE PRECISION   ALTSGN, ESTOLD, TEMP
  143: *     ..
  144: *     .. External Functions ..
  145:       INTEGER            IDAMAX
  146:       DOUBLE PRECISION   DASUM
  147:       EXTERNAL           IDAMAX, DASUM
  148: *     ..
  149: *     .. External Subroutines ..
  150:       EXTERNAL           DCOPY
  151: *     ..
  152: *     .. Intrinsic Functions ..
  153:       INTRINSIC          ABS, DBLE, NINT, SIGN
  154: *     ..
  155: *     .. Save statement ..
  156:       SAVE
  157: *     ..
  158: *     .. Executable Statements ..
  159: *
  160:       IF( KASE.EQ.0 ) THEN
  161:          DO 10 I = 1, N
  162:             X( I ) = ONE / DBLE( N )
  163:    10    CONTINUE
  164:          KASE = 1
  165:          JUMP = 1
  166:          RETURN
  167:       END IF
  168: *
  169:       GO TO ( 20, 40, 70, 110, 140 )JUMP
  170: *
  171: *     ................ ENTRY   (JUMP = 1)
  172: *     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY A*X.
  173: *
  174:    20 CONTINUE
  175:       IF( N.EQ.1 ) THEN
  176:          V( 1 ) = X( 1 )
  177:          EST = ABS( V( 1 ) )
  178: *        ... QUIT
  179:          GO TO 150
  180:       END IF
  181:       EST = DASUM( N, X, 1 )
  182: *
  183:       DO 30 I = 1, N
  184:          X( I ) = SIGN( ONE, X( I ) )
  185:          ISGN( I ) = NINT( X( I ) )
  186:    30 CONTINUE
  187:       KASE = 2
  188:       JUMP = 2
  189:       RETURN
  190: *
  191: *     ................ ENTRY   (JUMP = 2)
  192: *     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY TRANSPOSE(A)*X.
  193: *
  194:    40 CONTINUE
  195:       J = IDAMAX( N, X, 1 )
  196:       ITER = 2
  197: *
  198: *     MAIN LOOP - ITERATIONS 2,3,...,ITMAX.
  199: *
  200:    50 CONTINUE
  201:       DO 60 I = 1, N
  202:          X( I ) = ZERO
  203:    60 CONTINUE
  204:       X( J ) = ONE
  205:       KASE = 1
  206:       JUMP = 3
  207:       RETURN
  208: *
  209: *     ................ ENTRY   (JUMP = 3)
  210: *     X HAS BEEN OVERWRITTEN BY A*X.
  211: *
  212:    70 CONTINUE
  213:       CALL DCOPY( N, X, 1, V, 1 )
  214:       ESTOLD = EST
  215:       EST = DASUM( N, V, 1 )
  216:       DO 80 I = 1, N
  217:          IF( NINT( SIGN( ONE, X( I ) ) ).NE.ISGN( I ) )
  218:      $      GO TO 90
  219:    80 CONTINUE
  220: *     REPEATED SIGN VECTOR DETECTED, HENCE ALGORITHM HAS CONVERGED.
  221:       GO TO 120
  222: *
  223:    90 CONTINUE
  224: *     TEST FOR CYCLING.
  225:       IF( EST.LE.ESTOLD )
  226:      $   GO TO 120
  227: *
  228:       DO 100 I = 1, N
  229:          X( I ) = SIGN( ONE, X( I ) )
  230:          ISGN( I ) = NINT( X( I ) )
  231:   100 CONTINUE
  232:       KASE = 2
  233:       JUMP = 4
  234:       RETURN
  235: *
  236: *     ................ ENTRY   (JUMP = 4)
  237: *     X HAS BEEN OVERWRITTEN BY TRANSPOSE(A)*X.
  238: *
  239:   110 CONTINUE
  240:       JLAST = J
  241:       J = IDAMAX( N, X, 1 )
  242:       IF( ( X( JLAST ).NE.ABS( X( J ) ) ) .AND. ( ITER.LT.ITMAX ) ) THEN
  243:          ITER = ITER + 1
  244:          GO TO 50
  245:       END IF
  246: *
  247: *     ITERATION COMPLETE.  FINAL STAGE.
  248: *
  249:   120 CONTINUE
  250:       ALTSGN = ONE
  251:       DO 130 I = 1, N
  252:          X( I ) = ALTSGN*( ONE+DBLE( I-1 ) / DBLE( N-1 ) )
  253:          ALTSGN = -ALTSGN
  254:   130 CONTINUE
  255:       KASE = 1
  256:       JUMP = 5
  257:       RETURN
  258: *
  259: *     ................ ENTRY   (JUMP = 5)
  260: *     X HAS BEEN OVERWRITTEN BY A*X.
  261: *
  262:   140 CONTINUE
  263:       TEMP = TWO*( DASUM( N, X, 1 ) / DBLE( 3*N ) )
  264:       IF( TEMP.GT.EST ) THEN
  265:          CALL DCOPY( N, X, 1, V, 1 )
  266:          EST = TEMP
  267:       END IF
  268: *
  269:   150 CONTINUE
  270:       KASE = 0
  271:       RETURN
  272: *
  273: *     End of DLACON
  274: *
  275:       END

CVSweb interface <joel.bertrand@systella.fr>