Merge pull request #61 from cconrads-scicomp/improve-xSYEQUB

Improve xSYEQUB/xHEEQUB

6 files changed, 569 insertions, 579 deletions

diff --git a/SRC/cheequb.f b/SRC/cheequb.f
index 7f8dc20a..59b32c05 100644
--- a/SRC/cheequb.f
+++ b/SRC/cheequb.f
@@ -98,7 +98,7 @@
 *>
 *> \param[out] WORK
 *> \verbatim
-*>          WORK is COMPLEX array, dimension (3*N)
+*>          WORK is COMPLEX array, dimension (2*N)
 *> \endverbatim
 *>
 *> \param[out] INFO
@@ -151,14 +151,14 @@
 *
 *     .. Parameters ..
       REAL               ONE, ZERO
-      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
+      PARAMETER          ( ONE = 1.0E0, ZERO = 0.0E0 )
       INTEGER            MAX_ITER
       PARAMETER          ( MAX_ITER = 100 )
 *     ..
 *     .. Local Scalars ..
       INTEGER            I, J, ITER
-      REAL               AVG, STD, TOL, C0, C1, C2, T, U, SI, D,
-     $                   BASE, SMIN, SMAX, SMLNUM, BIGNUM, SCALE, SUMSQ
+      REAL               AVG, STD, TOL, C0, C1, C2, T, U, SI, D, BASE,
+     $                   SMIN, SMAX, SMLNUM, BIGNUM, SCALE, SUMSQ
       LOGICAL            UP
       COMPLEX            ZDUM
 *     ..
@@ -178,20 +178,22 @@
 *     ..
 *     .. Statement Function Definitions ..
       CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
+*     ..
+*     .. Executable Statements ..
 *
-*     Test input parameters.
+*     Test the input parameters.
 *
       INFO = 0
-      IF (.NOT. ( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) ) THEN
-        INFO = -1
+      IF ( .NOT. ( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) ) THEN
+         INFO = -1
       ELSE IF ( N .LT. 0 ) THEN
-        INFO = -2
+         INFO = -2
       ELSE IF ( LDA .LT. MAX( 1, N ) ) THEN
-        INFO = -4
+         INFO = -4
       END IF
       IF ( INFO .NE. 0 ) THEN
-        CALL XERBLA( 'CHEEQUB', -INFO )
-        RETURN
+         CALL XERBLA( 'CHEEQUB', -INFO )
+         RETURN
       END IF
 
       UP = LSAME( UPLO, 'U' )
@@ -200,12 +202,12 @@
 *     Quick return if possible.
 *
       IF ( N .EQ. 0 ) THEN
-        SCOND = ONE
-        RETURN
+         SCOND = ONE
+         RETURN
       END IF
 
       DO I = 1, N
-        S( I ) = ZERO
+         S( I ) = ZERO
       END DO
 
       AMAX = ZERO
@@ -226,102 +228,100 @@
             DO I = J+1, N
                S( I ) = MAX( S( I ), CABS1( A( I, J ) ) )
                S( J ) = MAX( S( J ), CABS1( A( I, J ) ) )
-               AMAX = MAX( AMAX, CABS1( A(I, J ) ) )
+               AMAX = MAX( AMAX, CABS1( A( I, J ) ) )
             END DO
          END DO
       END IF
       DO J = 1, N
-         S( J ) = 1.0 / S( J )
+         S( J ) = 1.0E0 / S( J )
       END DO
 
       TOL = ONE / SQRT( 2.0E0 * N )
 
       DO ITER = 1, MAX_ITER
-         SCALE = 0.0
-         SUMSQ = 0.0
-*       beta = |A|s
-        DO I = 1, N
-           WORK( I ) = ZERO
-        END DO
-        IF ( UP ) THEN
-           DO J = 1, N
-              DO I = 1, J-1
-                 T = CABS1( A( I, J ) )
-                 WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
-              END DO
-              WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
-           END DO
-        ELSE
-           DO J = 1, N
-              WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
-              DO I = J+1, N
-                 T = CABS1( A( I, J ) )
-                 WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
-              END DO
-           END DO
-        END IF
+         SCALE = 0.0E0
+         SUMSQ = 0.0E0
+*        beta = |A|s
+         DO I = 1, N
+            WORK( I ) = ZERO
+         END DO
+         IF ( UP ) THEN
+            DO J = 1, N
+               DO I = 1, J-1
+                  WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
+               END DO
+               WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
+            END DO
+         ELSE
+            DO J = 1, N
+               WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
+               DO I = J+1, N
+                  WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
+               END DO
+            END DO
+         END IF
 
-*       avg = s^T beta / n
-        AVG = 0.0
-        DO I = 1, N
-          AVG = AVG + S( I )*WORK( I )
-        END DO
-        AVG = AVG / N
+*        avg = s^T beta / n
+         AVG = 0.0E0
+         DO I = 1, N
+            AVG = AVG + S( I )*WORK( I )
+         END DO
+         AVG = AVG / N
 
-        STD = 0.0
-        DO I = 2*N+1, 3*N
-           WORK( I ) = S( I-2*N ) * WORK( I-2*N ) - AVG
-        END DO
-        CALL CLASSQ( N, WORK( 2*N+1 ), 1, SCALE, SUMSQ )
-        STD = SCALE * SQRT( SUMSQ / N )
+         STD = 0.0E0
+         DO I = N+1, 2*N
+            WORK( I ) = S( I-N ) * WORK( I-N ) - AVG
+         END DO
+         CALL CLASSQ( N, WORK( N+1 ), 1, SCALE, SUMSQ )
+         STD = SCALE * SQRT( SUMSQ / N )
 
-        IF ( STD .LT. TOL * AVG ) GOTO 999
+         IF ( STD .LT. TOL * AVG ) GOTO 999
 
-        DO I = 1, N
-          T = CABS1( A( I, I ) )
-          SI = S( I )
-          C2 = ( N-1 ) * T
-          C1 = ( N-2 ) * ( WORK( I ) - T*SI )
-          C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
+         DO I = 1, N
+            T = CABS1( A( I, I ) )
+            SI = S( I )
+            C2 = ( N-1 ) * T
+            C1 = ( N-2 ) * ( WORK( I ) - T*SI )
+            C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
+            D = C1*C1 - 4*C0*C2
 
-          D = C1*C1 - 4*C0*C2
-          IF ( D .LE. 0 ) THEN
-            INFO = -1
-            RETURN
-          END IF
-          SI = -2*C0 / ( C1 + SQRT( D ) )
+            IF ( D .LE. 0 ) THEN
+               INFO = -1
+               RETURN
+            END IF
+            SI = -2*C0 / ( C1 + SQRT( D ) )
 
-          D = SI - S(I)
-          U = ZERO
-          IF ( UP ) THEN
-            DO J = 1, I
-              T = CABS1( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = CABS1( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          ELSE
-            DO J = 1, I
-              T = CABS1( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = CABS1( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          END IF
-          AVG = AVG + ( U + WORK( I ) ) * D / N
-          S( I ) = SI
-        END DO
+            D = SI - S( I )
+            U = ZERO
+            IF ( UP ) THEN
+               DO J = 1, I
+                  T = CABS1( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = CABS1( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            ELSE
+               DO J = 1, I
+                  T = CABS1( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = CABS1( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            END IF
 
+            AVG = AVG + ( U + WORK( I ) ) * D / N
+            S( I ) = SI
+         END DO
       END DO
 
  999  CONTINUE
@@ -334,10 +334,10 @@
       BASE = SLAMCH( 'B' )
       U = ONE / LOG( BASE )
       DO I = 1, N
-        S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
-        SMIN = MIN( SMIN, S( I ) )
-        SMAX = MAX( SMAX, S( I ) )
+         S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
+         SMIN = MIN( SMIN, S( I ) )
+         SMAX = MAX( SMAX, S( I ) )
       END DO
       SCOND = MAX( SMIN, SMLNUM ) / MIN( SMAX, BIGNUM )
-
+*
       END
diff --git a/SRC/csyequb.f b/SRC/csyequb.f
index 5d54b5e8..7b6de6f1 100644
--- a/SRC/csyequb.f
+++ b/SRC/csyequb.f
@@ -98,7 +98,7 @@
 *>
 *> \param[out] WORK
 *> \verbatim
-*>          WORK is COMPLEX array, dimension (3*N)
+*>          WORK is COMPLEX array, dimension (2*N)
 *> \endverbatim
 *>
 *> \param[out] INFO
@@ -176,7 +176,7 @@
 *     .. Statement Functions ..
       REAL               CABS1
 *     ..
-*     Statement Function Definitions
+*     .. Statement Function Definitions ..
       CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) )
 *     ..
 *     .. Executable Statements ..
@@ -185,15 +185,15 @@
 *
       INFO = 0
       IF ( .NOT. ( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) ) THEN
-        INFO = -1
+         INFO = -1
       ELSE IF ( N .LT. 0 ) THEN
-        INFO = -2
+         INFO = -2
       ELSE IF ( LDA .LT. MAX( 1, N ) ) THEN
-        INFO = -4
+         INFO = -4
       END IF
       IF ( INFO .NE. 0 ) THEN
-        CALL XERBLA( 'CSYEQUB', -INFO )
-        RETURN
+         CALL XERBLA( 'CSYEQUB', -INFO )
+         RETURN
       END IF
 
       UP = LSAME( UPLO, 'U' )
@@ -202,12 +202,12 @@
 *     Quick return if possible.
 *
       IF ( N .EQ. 0 ) THEN
-        SCOND = ONE
-        RETURN
+         SCOND = ONE
+         RETURN
       END IF
 
       DO I = 1, N
-        S( I ) = ZERO
+         S( I ) = ZERO
       END DO
 
       AMAX = ZERO
@@ -218,7 +218,7 @@
                S( J ) = MAX( S( J ), CABS1( A( I, J ) ) )
                AMAX = MAX( AMAX, CABS1( A( I, J ) ) )
             END DO
-            S( J ) = MAX( S( J ), CABS1( A( J, J) ) )
+            S( J ) = MAX( S( J ), CABS1( A( J, J ) ) )
             AMAX = MAX( AMAX, CABS1( A( J, J ) ) )
          END DO
       ELSE
@@ -227,7 +227,7 @@
             AMAX = MAX( AMAX, CABS1( A( J, J ) ) )
             DO I = J+1, N
                S( I ) = MAX( S( I ), CABS1( A( I, J ) ) )
-               S( J ) = MAX( S( J ), CABS1 (A( I, J ) ) )
+               S( J ) = MAX( S( J ), CABS1( A( I, J ) ) )
                AMAX = MAX( AMAX, CABS1( A( I, J ) ) )
             END DO
          END DO
@@ -239,90 +239,89 @@
       TOL = ONE / SQRT( 2.0E0 * N )
 
       DO ITER = 1, MAX_ITER
-         SCALE = 0.0
-         SUMSQ = 0.0
-*       beta = |A|s
-        DO I = 1, N
-           WORK( I ) = ZERO
-        END DO
-        IF ( UP ) THEN
-           DO J = 1, N
-              DO I = 1, J-1
-                 T = CABS1( A( I, J ) )
-                 WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
-              END DO
-              WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
-           END DO
-        ELSE
-           DO J = 1, N
-              WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
-              DO I = J+1, N
-                 T = CABS1( A( I, J ) )
-                 WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
-              END DO
-           END DO
-        END IF
+         SCALE = 0.0E0
+         SUMSQ = 0.0E0
+*        beta = |A|s
+         DO I = 1, N
+             WORK( I ) = ZERO
+         END DO
+         IF ( UP ) THEN
+            DO J = 1, N
+               DO I = 1, J-1
+                  WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
+               END DO
+               WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
+            END DO
+         ELSE
+            DO J = 1, N
+               WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
+               DO I = J+1, N
+                  WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
+               END DO
+            END DO
+         END IF
 
-*       avg = s^T beta / n
-        AVG = 0.0
-        DO I = 1, N
-          AVG = AVG + S( I )*WORK( I )
-        END DO
-        AVG = AVG / N
+*        avg = s^T beta / n
+         AVG = 0.0E0
+         DO I = 1, N
+            AVG = AVG + S( I )*WORK( I )
+         END DO
+         AVG = AVG / N
 
-        STD = 0.0
-        DO I = N+1, 2*N
-           WORK( I ) = S( I-N ) * WORK( I-N ) - AVG
-        END DO
-        CALL CLASSQ( N, WORK( N+1 ), 1, SCALE, SUMSQ )
-        STD = SCALE * SQRT( SUMSQ / N )
+         STD = 0.0E0
+         DO I = N+1, 2*N
+            WORK( I ) = S( I-N ) * WORK( I-N ) - AVG
+         END DO
+         CALL CLASSQ( N, WORK( N+1 ), 1, SCALE, SUMSQ )
+         STD = SCALE * SQRT( SUMSQ / N )
 
-        IF ( STD .LT. TOL * AVG ) GOTO 999
+         IF ( STD .LT. TOL * AVG ) GOTO 999
 
-        DO I = 1, N
-          T = CABS1( A( I, I ) )
-          SI = S( I )
-          C2 = ( N-1 ) * T
-          C1 = ( N-2 ) * ( WORK( I ) - T*SI )
-          C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
-          D = C1*C1 - 4*C0*C2
+         DO I = 1, N
+            T = CABS1( A( I, I ) )
+            SI = S( I )
+            C2 = ( N-1 ) * T
+            C1 = ( N-2 ) * ( WORK( I ) - T*SI )
+            C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
+            D = C1*C1 - 4*C0*C2
 
-          IF ( D .LE. 0 ) THEN
-            INFO = -1
-            RETURN
-          END IF
-          SI = -2*C0 / ( C1 + SQRT( D ) )
+            IF ( D .LE. 0 ) THEN
+               INFO = -1
+               RETURN
+            END IF
+            SI = -2*C0 / ( C1 + SQRT( D ) )
 
-          D = SI - S( I )
-          U = ZERO
-          IF ( UP ) THEN
-            DO J = 1, I
-              T = CABS1( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = CABS1( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          ELSE
-            DO J = 1, I
-              T = CABS1( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = CABS1( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          END IF
-          AVG = AVG + ( U + WORK( I ) ) * D / N
-          S( I ) = SI
-        END DO
+            D = SI - S( I )
+            U = ZERO
+            IF ( UP ) THEN
+               DO J = 1, I
+                  T = CABS1( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = CABS1( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            ELSE
+               DO J = 1, I
+                  T = CABS1( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = CABS1( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            END IF
+
+            AVG = AVG + ( U + WORK( I ) ) * D / N
+            S( I ) = SI
+         END DO
       END DO
 
  999  CONTINUE
@@ -335,9 +334,9 @@
       BASE = SLAMCH( 'B' )
       U = ONE / LOG( BASE )
       DO I = 1, N
-        S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
-        SMIN = MIN( SMIN, S( I ) )
-        SMAX = MAX( SMAX, S( I ) )
+         S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
+         SMIN = MIN( SMIN, S( I ) )
+         SMAX = MAX( SMAX, S( I ) )
       END DO
       SCOND = MAX( SMIN, SMLNUM ) / MIN( SMAX, BIGNUM )
 *
diff --git a/SRC/dsyequb.f b/SRC/dsyequb.f
index 3c15cad4..7a549503 100644
--- a/SRC/dsyequb.f
+++ b/SRC/dsyequb.f
@@ -97,7 +97,7 @@
 *>
 *> \param[out] WORK
 *> \verbatim
-*>          WORK is DOUBLE PRECISION array, dimension (3*N)
+*>          WORK is DOUBLE PRECISION array, dimension (2*N)
 *> \endverbatim
 *>
 *> \param[out] INFO
@@ -149,7 +149,7 @@
 *
 *     .. Parameters ..
       DOUBLE PRECISION   ONE, ZERO
-      PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
+      PARAMETER          ( ONE = 1.0D0, ZERO = 0.0D0 )
       INTEGER            MAX_ITER
       PARAMETER          ( MAX_ITER = 100 )
 *     ..
@@ -172,19 +172,19 @@
 *     ..
 *     .. Executable Statements ..
 *
-*     Test input parameters.
+*     Test the input parameters.
 *
       INFO = 0
       IF ( .NOT. ( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) ) THEN
-        INFO = -1
+         INFO = -1
       ELSE IF ( N .LT. 0 ) THEN
-        INFO = -2
+         INFO = -2
       ELSE IF ( LDA .LT. MAX( 1, N ) ) THEN
-        INFO = -4
+         INFO = -4
       END IF
       IF ( INFO .NE. 0 ) THEN
-        CALL XERBLA( 'DSYEQUB', -INFO )
-        RETURN
+         CALL XERBLA( 'DSYEQUB', -INFO )
+         RETURN
       END IF
 
       UP = LSAME( UPLO, 'U' )
@@ -193,12 +193,12 @@
 *     Quick return if possible.
 *
       IF ( N .EQ. 0 ) THEN
-        SCOND = ONE
-        RETURN
+         SCOND = ONE
+         RETURN
       END IF
 
       DO I = 1, N
-        S( I ) = ZERO
+         S( I ) = ZERO
       END DO
 
       AMAX = ZERO
@@ -207,7 +207,7 @@
             DO I = 1, J-1
                S( I ) = MAX( S( I ), ABS( A( I, J ) ) )
                S( J ) = MAX( S( J ), ABS( A( I, J ) ) )
-               AMAX = MAX( AMAX, ABS( A(I, J) ) )
+               AMAX = MAX( AMAX, ABS( A( I, J ) ) )
             END DO
             S( J ) = MAX( S( J ), ABS( A( J, J ) ) )
             AMAX = MAX( AMAX, ABS( A( J, J ) ) )
@@ -224,99 +224,95 @@
          END DO
       END IF
       DO J = 1, N
-         S( J ) = 1.0D+0 / S( J )
+         S( J ) = 1.0D0 / S( J )
       END DO
 
-      TOL = ONE / SQRT(2.0D0 * N)
+      TOL = ONE / SQRT( 2.0D0 * N )
 
       DO ITER = 1, MAX_ITER
-         SCALE = 0.0D+0
-         SUMSQ = 0.0D+0
-*       BETA = |A|S
-        DO I = 1, N
-           WORK(I) = ZERO
-        END DO
-        IF ( UP ) THEN
-           DO J = 1, N
-              DO I = 1, J-1
-                 T = ABS( A( I, J ) )
-                 WORK( I ) = WORK( I ) + ABS( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + ABS( A( I, J ) ) * S( I )
-              END DO
-              WORK( J ) = WORK( J ) + ABS( A( J, J ) ) * S( J )
-           END DO
-        ELSE
-           DO J = 1, N
-              WORK( J ) = WORK( J ) + ABS( A( J, J ) ) * S( J )
-              DO I = J+1, N
-                 T = ABS( A( I, J ) )
-                 WORK( I ) = WORK( I ) + ABS( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + ABS( A( I, J ) ) * S( I )
-              END DO
-           END DO
-        END IF
-
-*       avg = s^T beta / n
-        AVG = 0.0D+0
-        DO I = 1, N
-          AVG = AVG + S( I )*WORK( I )
-        END DO
-        AVG = AVG / N
-
-        STD = 0.0D+0
-        DO I = 2*N+1, 3*N
-           WORK( I ) = S( I-2*N ) * WORK( I-2*N ) - AVG
-        END DO
-        CALL DLASSQ( N, WORK( 2*N+1 ), 1, SCALE, SUMSQ )
-        STD = SCALE * SQRT( SUMSQ / N )
+         SCALE = 0.0D0
+         SUMSQ = 0.0D0
+*        beta = |A|s
+         DO I = 1, N
+            WORK( I ) = ZERO
+         END DO
+         IF ( UP ) THEN
+            DO J = 1, N
+               DO I = 1, J-1
+                  WORK( I ) = WORK( I ) + ABS( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + ABS( A( I, J ) ) * S( I )
+               END DO
+               WORK( J ) = WORK( J ) + ABS( A( J, J ) ) * S( J )
+            END DO
+         ELSE
+            DO J = 1, N
+               WORK( J ) = WORK( J ) + ABS( A( J, J ) ) * S( J )
+               DO I = J+1, N
+                  WORK( I ) = WORK( I ) + ABS( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + ABS( A( I, J ) ) * S( I )
+               END DO
+            END DO
+         END IF
 
-        IF ( STD .LT. TOL * AVG ) GOTO 999
+*        avg = s^T beta / n
+         AVG = 0.0D0
+         DO I = 1, N
+            AVG = AVG + S( I )*WORK( I )
+         END DO
+         AVG = AVG / N
 
-        DO I = 1, N
-          T = ABS( A( I, I ) )
-          SI = S( I )
-          C2 = ( N-1 ) * T
-          C1 = ( N-2 ) * ( WORK( I ) - T*SI )
-          C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
-          D = C1*C1 - 4*C0*C2
+         STD = 0.0D0
+         DO I = N+1, 2*N
+            WORK( I ) = S( I-N ) * WORK( I-N ) - AVG
+         END DO
+         CALL DLASSQ( N, WORK( N+1 ), 1, SCALE, SUMSQ )
+         STD = SCALE * SQRT( SUMSQ / N )
 
-          IF ( D .LE. 0 ) THEN
-            INFO = -1
-            RETURN
-          END IF
-          SI = -2*C0 / ( C1 + SQRT( D ) )
+         IF ( STD .LT. TOL * AVG ) GOTO 999
 
-          D = SI - S( I )
-          U = ZERO
-          IF ( UP ) THEN
-            DO J = 1, I
-              T = ABS( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = ABS( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          ELSE
-            DO J = 1, I
-              T = ABS( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = ABS( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          END IF
+         DO I = 1, N
+            T = ABS( A( I, I ) )
+            SI = S( I )
+            C2 = ( N-1 ) * T
+            C1 = ( N-2 ) * ( WORK( I ) - T*SI )
+            C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
+            D = C1*C1 - 4*C0*C2
 
-          AVG = AVG + ( U + WORK( I ) ) * D / N
-          S( I ) = SI
+            IF ( D .LE. 0 ) THEN
+               INFO = -1
+               RETURN
+            END IF
+            SI = -2*C0 / ( C1 + SQRT( D ) )
 
-        END DO
+            D = SI - S( I )
+            U = ZERO
+            IF ( UP ) THEN
+               DO J = 1, I
+                  T = ABS( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = ABS( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            ELSE
+               DO J = 1, I
+                  T = ABS( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = ABS( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            END IF
 
+            AVG = AVG + ( U + WORK( I ) ) * D / N
+            S( I ) = SI
+         END DO
       END DO
 
  999  CONTINUE
@@ -325,13 +321,13 @@
       BIGNUM = ONE / SMLNUM
       SMIN = BIGNUM
       SMAX = ZERO
-      T = ONE / SQRT(AVG)
+      T = ONE / SQRT( AVG )
       BASE = DLAMCH( 'B' )
       U = ONE / LOG( BASE )
       DO I = 1, N
-        S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
-        SMIN = MIN( SMIN, S( I ) )
-        SMAX = MAX( SMAX, S( I ) )
+         S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
+         SMIN = MIN( SMIN, S( I ) )
+         SMAX = MAX( SMAX, S( I ) )
       END DO
       SCOND = MAX( SMIN, SMLNUM ) / MIN( SMAX, BIGNUM )
 *
diff --git a/SRC/ssyequb.f b/SRC/ssyequb.f
index eb143906..56dfca84 100644
--- a/SRC/ssyequb.f
+++ b/SRC/ssyequb.f
@@ -97,7 +97,7 @@
 *>
 *> \param[out] WORK
 *> \verbatim
-*>          WORK is REAL array, dimension (3*N)
+*>          WORK is REAL array, dimension (2*N)
 *> \endverbatim
 *>
 *> \param[out] INFO
@@ -149,7 +149,7 @@
 *
 *     .. Parameters ..
       REAL               ONE, ZERO
-      PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
+      PARAMETER          ( ONE = 1.0E0, ZERO = 0.0E0 )
       INTEGER            MAX_ITER
       PARAMETER          ( MAX_ITER = 100 )
 *     ..
@@ -172,19 +172,19 @@
 *     ..
 *     .. Executable Statements ..
 *
-*     Test input parameters.
+*     Test the input parameters.
 *
       INFO = 0
       IF ( .NOT. ( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) ) THEN
-        INFO = -1
+         INFO = -1
       ELSE IF ( N .LT. 0 ) THEN
-        INFO = -2
+         INFO = -2
       ELSE IF ( LDA .LT. MAX( 1, N ) ) THEN
-        INFO = -4
+         INFO = -4
       END IF
       IF ( INFO .NE. 0 ) THEN
-        CALL XERBLA( 'SSYEQUB', -INFO )
-        RETURN
+         CALL XERBLA( 'SSYEQUB', -INFO )
+         RETURN
       END IF
 
       UP = LSAME( UPLO, 'U' )
@@ -193,12 +193,12 @@
 *     Quick return if possible.
 *
       IF ( N .EQ. 0 ) THEN
-        SCOND = ONE
-        RETURN
+         SCOND = ONE
+         RETURN
       END IF
 
       DO I = 1, N
-        S( I ) = ZERO
+         S( I ) = ZERO
       END DO
 
       AMAX = ZERO
@@ -207,7 +207,7 @@
             DO I = 1, J-1
                S( I ) = MAX( S( I ), ABS( A( I, J ) ) )
                S( J ) = MAX( S( J ), ABS( A( I, J ) ) )
-               AMAX = MAX( AMAX, ABS( A(I, J) ) )
+               AMAX = MAX( AMAX, ABS( A( I, J ) ) )
             END DO
             S( J ) = MAX( S( J ), ABS( A( J, J ) ) )
             AMAX = MAX( AMAX, ABS( A( J, J ) ) )
@@ -224,99 +224,95 @@
          END DO
       END IF
       DO J = 1, N
-         S( J ) = 1.0 / S( J )
+         S( J ) = 1.0E0 / S( J )
       END DO
 
-      TOL = ONE / SQRT(2.0E0 * N)
+      TOL = ONE / SQRT( 2.0E0 * N )
 
       DO ITER = 1, MAX_ITER
-         SCALE = 0.0
-         SUMSQ = 0.0
-*       BETA = |A|S
-        DO I = 1, N
-           WORK(I) = ZERO
-        END DO
-        IF ( UP ) THEN
-           DO J = 1, N
-              DO I = 1, J-1
-                 T = ABS( A( I, J ) )
-                 WORK( I ) = WORK( I ) + ABS( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + ABS( A( I, J ) ) * S( I )
-              END DO
-              WORK( J ) = WORK( J ) + ABS( A( J, J ) ) * S( J )
-           END DO
-        ELSE
-           DO J = 1, N
-              WORK( J ) = WORK( J ) + ABS( A( J, J ) ) * S( J )
-              DO I = J+1, N
-                 T = ABS( A( I, J ) )
-                 WORK( I ) = WORK( I ) + ABS( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + ABS( A( I, J ) ) * S( I )
-              END DO
-           END DO
-        END IF
-
-*       avg = s^T beta / n
-        AVG = 0.0
-        DO I = 1, N
-          AVG = AVG + S( I )*WORK( I )
-        END DO
-        AVG = AVG / N
-
-        STD = 0.0
-        DO I = 2*N+1, 3*N
-           WORK( I ) = S( I-2*N ) * WORK( I-2*N ) - AVG
-        END DO
-        CALL SLASSQ( N, WORK( 2*N+1 ), 1, SCALE, SUMSQ )
-        STD = SCALE * SQRT( SUMSQ / N )
+         SCALE = 0.0E0
+         SUMSQ = 0.0E0
+*        beta = |A|s
+         DO I = 1, N
+            WORK( I ) = ZERO
+         END DO
+         IF ( UP ) THEN
+            DO J = 1, N
+               DO I = 1, J-1
+                  WORK( I ) = WORK( I ) + ABS( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + ABS( A( I, J ) ) * S( I )
+               END DO
+               WORK( J ) = WORK( J ) + ABS( A( J, J ) ) * S( J )
+            END DO
+         ELSE
+            DO J = 1, N
+               WORK( J ) = WORK( J ) + ABS( A( J, J ) ) * S( J )
+               DO I = J+1, N
+                  WORK( I ) = WORK( I ) + ABS( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + ABS( A( I, J ) ) * S( I )
+               END DO
+            END DO
+         END IF
 
-        IF ( STD .LT. TOL * AVG ) GOTO 999
+*        avg = s^T beta / n
+         AVG = 0.0E0
+         DO I = 1, N
+            AVG = AVG + S( I )*WORK( I )
+         END DO
+         AVG = AVG / N
 
-        DO I = 1, N
-          T = ABS( A( I, I ) )
-          SI = S( I )
-          C2 = ( N-1 ) * T
-          C1 = ( N-2 ) * ( WORK( I ) - T*SI )
-          C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
-          D = C1*C1 - 4*C0*C2
+         STD = 0.0E0
+         DO I = N+1, 2*N
+            WORK( I ) = S( I-N ) * WORK( I-N ) - AVG
+         END DO
+         CALL SLASSQ( N, WORK( N+1 ), 1, SCALE, SUMSQ )
+         STD = SCALE * SQRT( SUMSQ / N )
 
-          IF ( D .LE. 0 ) THEN
-            INFO = -1
-            RETURN
-          END IF
-          SI = -2*C0 / ( C1 + SQRT( D ) )
+         IF ( STD .LT. TOL * AVG ) GOTO 999
 
-          D = SI - S( I )
-          U = ZERO
-          IF ( UP ) THEN
-            DO J = 1, I
-              T = ABS( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = ABS( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          ELSE
-            DO J = 1, I
-              T = ABS( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = ABS( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          END IF
+         DO I = 1, N
+            T = ABS( A( I, I ) )
+            SI = S( I )
+            C2 = ( N-1 ) * T
+            C1 = ( N-2 ) * ( WORK( I ) - T*SI )
+            C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
+            D = C1*C1 - 4*C0*C2
 
-          AVG = AVG + ( U + WORK( I ) ) * D / N
-          S( I ) = SI
+            IF ( D .LE. 0 ) THEN
+               INFO = -1
+               RETURN
+            END IF
+            SI = -2*C0 / ( C1 + SQRT( D ) )
 
-        END DO
+            D = SI - S( I )
+            U = ZERO
+            IF ( UP ) THEN
+               DO J = 1, I
+                  T = ABS( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = ABS( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            ELSE
+               DO J = 1, I
+                  T = ABS( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = ABS( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            END IF
 
+            AVG = AVG + ( U + WORK( I ) ) * D / N
+            S( I ) = SI
+         END DO
       END DO
 
  999  CONTINUE
@@ -325,13 +321,13 @@
       BIGNUM = ONE / SMLNUM
       SMIN = BIGNUM
       SMAX = ZERO
-      T = ONE / SQRT(AVG)
+      T = ONE / SQRT( AVG )
       BASE = SLAMCH( 'B' )
       U = ONE / LOG( BASE )
       DO I = 1, N
-        S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
-        SMIN = MIN( SMIN, S( I ) )
-        SMAX = MAX( SMAX, S( I ) )
+         S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
+         SMIN = MIN( SMIN, S( I ) )
+         SMAX = MAX( SMAX, S( I ) )
       END DO
       SCOND = MAX( SMIN, SMLNUM ) / MIN( SMAX, BIGNUM )
 *
diff --git a/SRC/zheequb.f b/SRC/zheequb.f
index 341b7d22..1e06299a 100644
--- a/SRC/zheequb.f
+++ b/SRC/zheequb.f
@@ -98,7 +98,7 @@
 *>
 *> \param[out] WORK
 *> \verbatim
-*>          WORK is COMPLEX*16 array, dimension (3*N)
+*>          WORK is COMPLEX*16 array, dimension (2*N)
 *> \endverbatim
 *>
 *> \param[out] INFO
@@ -151,14 +151,14 @@
 *
 *     .. Parameters ..
       DOUBLE PRECISION   ONE, ZERO
-      PARAMETER          ( ONE = 1.0D+0, ZERO = 0.0D+0 )
+      PARAMETER          ( ONE = 1.0D0, ZERO = 0.0D0 )
       INTEGER            MAX_ITER
       PARAMETER          ( MAX_ITER = 100 )
 *     ..
 *     .. Local Scalars ..
       INTEGER            I, J, ITER
-      DOUBLE PRECISION   AVG, STD, TOL, C0, C1, C2, T, U, SI, D,
-     $                   BASE, SMIN, SMAX, SMLNUM, BIGNUM, SCALE, SUMSQ
+      DOUBLE PRECISION   AVG, STD, TOL, C0, C1, C2, T, U, SI, D, BASE,
+     $                   SMIN, SMAX, SMLNUM, BIGNUM, SCALE, SUMSQ
       LOGICAL            UP
       COMPLEX*16         ZDUM
 *     ..
@@ -178,20 +178,22 @@
 *     ..
 *     .. Statement Function Definitions ..
       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
+*     ..
+*     .. Executable Statements ..
 *
-*     Test input parameters.
+*     Test the input parameters.
 *
       INFO = 0
-      IF (.NOT. ( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) ) THEN
-        INFO = -1
+      IF ( .NOT. ( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) ) THEN
+         INFO = -1
       ELSE IF ( N .LT. 0 ) THEN
-        INFO = -2
+         INFO = -2
       ELSE IF ( LDA .LT. MAX( 1, N ) ) THEN
-        INFO = -4
+         INFO = -4
       END IF
       IF ( INFO .NE. 0 ) THEN
-        CALL XERBLA( 'ZHEEQUB', -INFO )
-        RETURN
+         CALL XERBLA( 'ZHEEQUB', -INFO )
+         RETURN
       END IF
 
       UP = LSAME( UPLO, 'U' )
@@ -200,12 +202,12 @@
 *     Quick return if possible.
 *
       IF ( N .EQ. 0 ) THEN
-        SCOND = ONE
-        RETURN
+         SCOND = ONE
+         RETURN
       END IF
 
       DO I = 1, N
-        S( I ) = ZERO
+         S( I ) = ZERO
       END DO
 
       AMAX = ZERO
@@ -226,102 +228,100 @@
             DO I = J+1, N
                S( I ) = MAX( S( I ), CABS1( A( I, J ) ) )
                S( J ) = MAX( S( J ), CABS1( A( I, J ) ) )
-               AMAX = MAX( AMAX, CABS1( A(I, J ) ) )
+               AMAX = MAX( AMAX, CABS1( A( I, J ) ) )
             END DO
          END DO
       END IF
       DO J = 1, N
-         S( J ) = 1.0D+0 / S( J )
+         S( J ) = 1.0D0 / S( J )
       END DO
 
       TOL = ONE / SQRT( 2.0D0 * N )
 
       DO ITER = 1, MAX_ITER
-         SCALE = 0.0D+0
-         SUMSQ = 0.0D+0
-*       beta = |A|s
-        DO I = 1, N
-           WORK( I ) = ZERO
-        END DO
-        IF ( UP ) THEN
-           DO J = 1, N
-              DO I = 1, J-1
-                 T = CABS1( A( I, J ) )
-                 WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
-              END DO
-              WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
-           END DO
-        ELSE
-           DO J = 1, N
-              WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
-              DO I = J+1, N
-                 T = CABS1( A( I, J ) )
-                 WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
-              END DO
-           END DO
-        END IF
+         SCALE = 0.0D0
+         SUMSQ = 0.0D0
+*        beta = |A|s
+         DO I = 1, N
+            WORK( I ) = ZERO
+         END DO
+         IF ( UP ) THEN
+            DO J = 1, N
+               DO I = 1, J-1
+                  WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
+               END DO
+               WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
+            END DO
+         ELSE
+            DO J = 1, N
+               WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
+               DO I = J+1, N
+                  WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
+               END DO
+            END DO
+         END IF
 
-*       avg = s^T beta / n
-        AVG = 0.0D+0
-        DO I = 1, N
-          AVG = AVG + S( I )*WORK( I )
-        END DO
-        AVG = AVG / N
+*        avg = s^T beta / n
+         AVG = 0.0D0
+         DO I = 1, N
+            AVG = AVG + S( I )*WORK( I )
+         END DO
+         AVG = AVG / N
 
-        STD = 0.0D+0
-        DO I = 2*N+1, 3*N
-           WORK( I ) = S( I-2*N ) * WORK( I-2*N ) - AVG
-        END DO
-        CALL ZLASSQ( N, WORK( 2*N+1 ), 1, SCALE, SUMSQ )
-        STD = SCALE * SQRT( SUMSQ / N )
+         STD = 0.0D0
+         DO I = N+1, N
+            WORK( I ) = S( I-N ) * WORK( I-N ) - AVG
+         END DO
+         CALL ZLASSQ( N, WORK( N+1 ), 1, SCALE, SUMSQ )
+         STD = SCALE * SQRT( SUMSQ / N )
 
-        IF ( STD .LT. TOL * AVG ) GOTO 999
+         IF ( STD .LT. TOL * AVG ) GOTO 999
 
-        DO I = 1, N
-          T = CABS1( A( I, I ) )
-          SI = S( I )
-          C2 = ( N-1 ) * T
-          C1 = ( N-2 ) * ( WORK( I ) - T*SI )
-          C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
+         DO I = 1, N
+            T = CABS1( A( I, I ) )
+            SI = S( I )
+            C2 = ( N-1 ) * T
+            C1 = ( N-2 ) * ( WORK( I ) - T*SI )
+            C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
+            D = C1*C1 - 4*C0*C2
 
-          D = C1*C1 - 4*C0*C2
-          IF ( D .LE. 0 ) THEN
-            INFO = -1
-            RETURN
-          END IF
-          SI = -2*C0 / ( C1 + SQRT( D ) )
+            IF ( D .LE. 0 ) THEN
+               INFO = -1
+               RETURN
+            END IF
+            SI = -2*C0 / ( C1 + SQRT( D ) )
 
-          D = SI - S(I)
-          U = ZERO
-          IF ( UP ) THEN
-            DO J = 1, I
-              T = CABS1( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = CABS1( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          ELSE
-            DO J = 1, I
-              T = CABS1( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = CABS1( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          END IF
-          AVG = AVG + ( U + WORK( I ) ) * D / N
-          S( I ) = SI
-        END DO
+            D = SI - S( I )
+            U = ZERO
+            IF ( UP ) THEN
+               DO J = 1, I
+                  T = CABS1( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = CABS1( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            ELSE
+               DO J = 1, I
+                  T = CABS1( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = CABS1( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            END IF
 
+            AVG = AVG + ( U + WORK( I ) ) * D / N
+            S( I ) = SI
+         END DO
       END DO
 
  999  CONTINUE
@@ -334,10 +334,10 @@
       BASE = DLAMCH( 'B' )
       U = ONE / LOG( BASE )
       DO I = 1, N
-        S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
-        SMIN = MIN( SMIN, S( I ) )
-        SMAX = MAX( SMAX, S( I ) )
+         S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
+         SMIN = MIN( SMIN, S( I ) )
+         SMAX = MAX( SMAX, S( I ) )
       END DO
       SCOND = MAX( SMIN, SMLNUM ) / MIN( SMAX, BIGNUM )
-
+*
       END
diff --git a/SRC/zsyequb.f b/SRC/zsyequb.f
index 8514388e..847f5064 100644
--- a/SRC/zsyequb.f
+++ b/SRC/zsyequb.f
@@ -98,7 +98,7 @@
 *>
 *> \param[out] WORK
 *> \verbatim
-*>          WORK is COMPLEX*16 array, dimension (3*N)
+*>          WORK is COMPLEX*16 array, dimension (2*N)
 *> \endverbatim
 *>
 *> \param[out] INFO
@@ -176,7 +176,7 @@
 *     .. Statement Functions ..
       DOUBLE PRECISION   CABS1
 *     ..
-*     Statement Function Definitions
+*     .. Statement Function Definitions ..
       CABS1( ZDUM ) = ABS( DBLE( ZDUM ) ) + ABS( DIMAG( ZDUM ) )
 *     ..
 *     .. Executable Statements ..
@@ -185,15 +185,15 @@
 *
       INFO = 0
       IF ( .NOT. ( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) ) THEN
-        INFO = -1
+         INFO = -1
       ELSE IF ( N .LT. 0 ) THEN
-        INFO = -2
+         INFO = -2
       ELSE IF ( LDA .LT. MAX( 1, N ) ) THEN
-        INFO = -4
+         INFO = -4
       END IF
       IF ( INFO .NE. 0 ) THEN
-        CALL XERBLA( 'ZSYEQUB', -INFO )
-        RETURN
+         CALL XERBLA( 'ZSYEQUB', -INFO )
+         RETURN
       END IF
 
       UP = LSAME( UPLO, 'U' )
@@ -202,12 +202,12 @@
 *     Quick return if possible.
 *
       IF ( N .EQ. 0 ) THEN
-        SCOND = ONE
-        RETURN
+         SCOND = ONE
+         RETURN
       END IF
 
       DO I = 1, N
-        S( I ) = ZERO
+         S( I ) = ZERO
       END DO
 
       AMAX = ZERO
@@ -218,7 +218,7 @@
                S( J ) = MAX( S( J ), CABS1( A( I, J ) ) )
                AMAX = MAX( AMAX, CABS1( A( I, J ) ) )
             END DO
-            S( J ) = MAX( S( J ), CABS1( A( J, J) ) )
+            S( J ) = MAX( S( J ), CABS1( A( J, J ) ) )
             AMAX = MAX( AMAX, CABS1( A( J, J ) ) )
          END DO
       ELSE
@@ -227,102 +227,101 @@
             AMAX = MAX( AMAX, CABS1( A( J, J ) ) )
             DO I = J+1, N
                S( I ) = MAX( S( I ), CABS1( A( I, J ) ) )
-               S( J ) = MAX( S( J ), CABS1 (A( I, J ) ) )
+               S( J ) = MAX( S( J ), CABS1( A( I, J ) ) )
                AMAX = MAX( AMAX, CABS1( A( I, J ) ) )
             END DO
          END DO
       END IF
       DO J = 1, N
-         S( J ) = 1.0D+0 / S( J )
+         S( J ) = 1.0D0 / S( J )
       END DO
 
       TOL = ONE / SQRT( 2.0D0 * N )
 
       DO ITER = 1, MAX_ITER
-         SCALE = 0.0D+0
-         SUMSQ = 0.0D+0
-*       beta = |A|s
-        DO I = 1, N
-           WORK( I ) = ZERO
-        END DO
-        IF ( UP ) THEN
-           DO J = 1, N
-              DO I = 1, J-1
-                 T = CABS1( A( I, J ) )
-                 WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
-              END DO
-              WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
-           END DO
-        ELSE
-           DO J = 1, N
-              WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
-              DO I = J+1, N
-                 T = CABS1( A( I, J ) )
-                 WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
-                 WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
-              END DO
-           END DO
-        END IF
+         SCALE = 0.0D0
+         SUMSQ = 0.0D0
+*        beta = |A|s
+         DO I = 1, N
+            WORK( I ) = ZERO
+         END DO
+         IF ( UP ) THEN
+            DO J = 1, N
+               DO I = 1, J-1
+                  WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
+               END DO
+               WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
+            END DO
+         ELSE
+            DO J = 1, N
+               WORK( J ) = WORK( J ) + CABS1( A( J, J ) ) * S( J )
+               DO I = J+1, N
+                  WORK( I ) = WORK( I ) + CABS1( A( I, J ) ) * S( J )
+                  WORK( J ) = WORK( J ) + CABS1( A( I, J ) ) * S( I )
+               END DO
+            END DO
+         END IF
 
-*       avg = s^T beta / n
-        AVG = 0.0D+0
-        DO I = 1, N
-          AVG = AVG + S( I )*WORK( I )
-        END DO
-        AVG = AVG / N
+*        avg = s^T beta / n
+         AVG = 0.0D0
+         DO I = 1, N
+            AVG = AVG + S( I )*WORK( I )
+         END DO
+         AVG = AVG / N
 
-        STD = 0.0D+0
-        DO I = N+1, 2*N
-           WORK( I ) = S( I-N ) * WORK( I-N ) - AVG
-        END DO
-        CALL ZLASSQ( N, WORK( N+1 ), 1, SCALE, SUMSQ )
-        STD = SCALE * SQRT( SUMSQ / N )
+         STD = 0.0D0
+         DO I = N+1, 2*N
+            WORK( I ) = S( I-N ) * WORK( I-N ) - AVG
+         END DO
+         CALL ZLASSQ( N, WORK( N+1 ), 1, SCALE, SUMSQ )
+         STD = SCALE * SQRT( SUMSQ / N )
 
-        IF ( STD .LT. TOL * AVG ) GOTO 999
+         IF ( STD .LT. TOL * AVG ) GOTO 999
 
-        DO I = 1, N
-          T = CABS1( A( I, I ) )
-          SI = S( I )
-          C2 = ( N-1 ) * T
-          C1 = ( N-2 ) * ( WORK( I ) - T*SI )
-          C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
-          D = C1*C1 - 4*C0*C2
+         DO I = 1, N
+            T = CABS1( A( I, I ) )
+            SI = S( I )
+            C2 = ( N-1 ) * T
+            C1 = ( N-2 ) * ( WORK( I ) - T*SI )
+            C0 = -(T*SI)*SI + 2*WORK( I )*SI - N*AVG
+            D = C1*C1 - 4*C0*C2
 
-          IF ( D .LE. 0 ) THEN
-            INFO = -1
-            RETURN
-          END IF
-          SI = -2*C0 / ( C1 + SQRT( D ) )
+            IF ( D .LE. 0 ) THEN
+               INFO = -1
+               RETURN
+            END IF
+            SI = -2*C0 / ( C1 + SQRT( D ) )
 
-          D = SI - S( I )
-          U = ZERO
-          IF ( UP ) THEN
-            DO J = 1, I
-              T = CABS1( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = CABS1( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          ELSE
-            DO J = 1, I
-              T = CABS1( A( I, J ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-            DO J = I+1,N
-              T = CABS1( A( J, I ) )
-              U = U + S( J )*T
-              WORK( J ) = WORK( J ) + D*T
-            END DO
-          END IF
-          AVG = AVG + ( U + WORK( I ) ) * D / N
-          S( I ) = SI
-        END DO
+            D = SI - S( I )
+            U = ZERO
+            IF ( UP ) THEN
+               DO J = 1, I
+                  T = CABS1( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = CABS1( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            ELSE
+               DO J = 1, I
+                  T = CABS1( A( I, J ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+               DO J = I+1,N
+                  T = CABS1( A( J, I ) )
+                  U = U + S( J )*T
+                  WORK( J ) = WORK( J ) + D*T
+               END DO
+            END IF
+
+            AVG = AVG + ( U + WORK( I ) ) * D / N
+            S( I ) = SI
+         END DO
       END DO
 
  999  CONTINUE
@@ -335,9 +334,9 @@
       BASE = DLAMCH( 'B' )
       U = ONE / LOG( BASE )
       DO I = 1, N
-        S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
-        SMIN = MIN( SMIN, S( I ) )
-        SMAX = MAX( SMAX, S( I ) )
+         S( I ) = BASE ** INT( U * LOG( S( I ) * T ) )
+         SMIN = MIN( SMIN, S( I ) )
+         SMAX = MAX( SMAX, S( I ) )
       END DO
       SCOND = MAX( SMIN, SMLNUM ) / MIN( SMAX, BIGNUM )
 *