slacn2()

subroutine slacn2	(	integer	n,
		real, dimension( * )	v,
		real, dimension( * )	x,
		integer, dimension( * )	isgn,
		real	est,
		integer	kase,
		integer, dimension( 3 )	isave
	)

SLACN2 estimates the 1-norm of a square matrix, using reverse communication for evaluating matrix-vector products.

Download SLACN2 + dependencies [TGZ] [ZIP] [TXT]

Purpose:

 SLACN2 estimates the 1-norm of a square, real matrix A.
 Reverse communication is used for evaluating matrix-vector products.

Parameters

[in]	N	N is INTEGER The order of the matrix. N >= 1.
[out]	V	V is REAL array, dimension (N) On the final return, V = A*W, where EST = norm(V)/norm(W) (W is not returned).
[in,out]	X	X is REAL array, dimension (N) On an intermediate return, X should be overwritten by A * X, if KASE=1, A**T * X, if KASE=2, and SLACN2 must be re-called with all the other parameters unchanged.
[out]	ISGN	ISGN is INTEGER array, dimension (N)
[in,out]	EST	EST is REAL On entry with KASE = 1 or 2 and ISAVE(1) = 3, EST should be unchanged from the previous call to SLACN2. On exit, EST is an estimate (a lower bound) for norm(A).
[in,out]	KASE	KASE is INTEGER On the initial call to SLACN2, KASE should be 0. On an intermediate return, KASE will be 1 or 2, indicating whether X should be overwritten by A * X or A**T * X. On the final return from SLACN2, KASE will again be 0.
[in,out]	ISAVE	ISAVE is INTEGER array, dimension (3) ISAVE is used to save variables between calls to SLACN2

Author

Univ. of Tennessee

Univ. of California Berkeley

Univ. of Colorado Denver

NAG Ltd.

Further Details:

  Originally named SONEST, dated March 16, 1988.

  This is a thread safe version of SLACON, which uses the array ISAVE
  in place of a SAVE statement, as follows:

     SLACON     SLACN2
      JUMP     ISAVE(1)
      J        ISAVE(2)
      ITER     ISAVE(3)

Contributors:

Nick Higham, University of Manchester

References:

N.J. Higham, "FORTRAN codes for estimating the one-norm of a real or complex matrix, with applications to condition estimation", ACM Trans. Math. Soft., vol. 14, no. 4, pp. 381-396, December 1988.

Definition at line 135 of file slacn2.f.

*
*  -- LAPACK auxiliary routine --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
*     .. Scalar Arguments ..
      INTEGER            KASE, N
      REAL               EST
*     ..
*     .. Array Arguments ..
      INTEGER            ISGN( * ), ISAVE( 3 )
      REAL               V( * ), X( * )
*     ..
*
*  =====================================================================
*
*     .. Parameters ..
      INTEGER            ITMAX
      parameter( itmax = 5 )
      REAL               ZERO, ONE, TWO
      parameter( zero = 0.0e+0, one = 1.0e+0, two = 2.0e+0 )
*     ..
*     .. Local Scalars ..
      INTEGER            I, JLAST
      REAL               ALTSGN, ESTOLD, TEMP, XS
*     ..
*     .. External Functions ..
      INTEGER            ISAMAX
      REAL               SASUM
      EXTERNAL           isamax, sasum
*     ..
*     .. External Subroutines ..
      EXTERNAL           scopy
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          abs, nint, real
*     ..
*     .. Executable Statements ..
*
      IF( kase.EQ.0 ) THEN
         DO 10 i = 1, n
            x( i ) = one / real( n )
   10    CONTINUE
         kase = 1
         isave( 1 ) = 1
         RETURN
      END IF
*
      GO TO ( 20, 40, 70, 110, 140 )isave( 1 )
*
*     ................ ENTRY   (ISAVE( 1 ) = 1)
*     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY A*X.
*
   20 CONTINUE
      IF( n.EQ.1 ) THEN
         v( 1 ) = x( 1 )
         est = abs( v( 1 ) )
*        ... QUIT
         GO TO 150
      END IF
      est = sasum( n, x, 1 )
*
      DO 30 i = 1, n
         IF( x(i).GE.zero ) THEN
            x(i) = one
         ELSE
            x(i) = -one
         END IF
         isgn( i ) = nint( x( i ) )
   30 CONTINUE
      kase = 2
      isave( 1 ) = 2
      RETURN
*
*     ................ ENTRY   (ISAVE( 1 ) = 2)
*     FIRST ITERATION.  X HAS BEEN OVERWRITTEN BY TRANSPOSE(A)*X.
*
   40 CONTINUE
      isave( 2 ) = isamax( n, x, 1 )
      isave( 3 ) = 2
*
*     MAIN LOOP - ITERATIONS 2,3,...,ITMAX.
*
   50 CONTINUE
      DO 60 i = 1, n
         x( i ) = zero
   60 CONTINUE
      x( isave( 2 ) ) = one
      kase = 1
      isave( 1 ) = 3
      RETURN
*
*     ................ ENTRY   (ISAVE( 1 ) = 3)
*     X HAS BEEN OVERWRITTEN BY A*X.
*
   70 CONTINUE
      CALL scopy( n, x, 1, v, 1 )
      estold = est
      est = sasum( n, v, 1 )
      DO 80 i = 1, n
         IF( x(i).GE.zero ) THEN
            xs = one
         ELSE
            xs = -one
         END IF
         IF( nint( xs ).NE.isgn( i ) )
     $      GO TO 90
   80 CONTINUE
*     REPEATED SIGN VECTOR DETECTED, HENCE ALGORITHM HAS CONVERGED.
      GO TO 120
*
   90 CONTINUE
*     TEST FOR CYCLING.
      IF( est.LE.estold )
     $   GO TO 120
*
      DO 100 i = 1, n
         IF( x(i).GE.zero ) THEN
            x(i) = one
         ELSE
            x(i) = -one
         END IF
         isgn( i ) = nint( x( i ) )
  100 CONTINUE
      kase = 2
      isave( 1 ) = 4
      RETURN
*
*     ................ ENTRY   (ISAVE( 1 ) = 4)
*     X HAS BEEN OVERWRITTEN BY TRANSPOSE(A)*X.
*
  110 CONTINUE
      jlast = isave( 2 )
      isave( 2 ) = isamax( n, x, 1 )
      IF( ( x( jlast ).NE.abs( x( isave( 2 ) ) ) ) .AND.
     $    ( isave( 3 ).LT.itmax ) ) THEN
         isave( 3 ) = isave( 3 ) + 1
         GO TO 50
      END IF
*
*     ITERATION COMPLETE.  FINAL STAGE.
*
  120 CONTINUE
      altsgn = one
      DO 130 i = 1, n
         x( i ) = altsgn*( one+real( i-1 ) / real( n-1 ) )
         altsgn = -altsgn
  130 CONTINUE
      kase = 1
      isave( 1 ) = 5
      RETURN
*
*     ................ ENTRY   (ISAVE( 1 ) = 5)
*     X HAS BEEN OVERWRITTEN BY A*X.
*
  140 CONTINUE
      temp = two*( sasum( n, x, 1 ) / real( 3*n ) )
      IF( temp.GT.est ) THEN
         CALL scopy( n, x, 1, v, 1 )
         est = temp
      END IF
*
  150 CONTINUE
      kase = 0
      RETURN
*
*     End of SLACN2
*

Here is the call graph for this function:

Here is the caller graph for this function: