Now using 2-norm to compute vector norms of row and column for balancing algorithm.

This seems to fix problems with balancing causing very large backward error for certain
Hessenberg matrices, including the often cited example of Watkins.

4 files changed, 21 insertions, 44 deletions

diff --git a/SRC/cgebal.f b/SRC/cgebal.f
index 4bb5a2f0..47207ea3 100644
--- a/SRC/cgebal.f
+++ b/SRC/cgebal.f
@@ -195,8 +195,8 @@
 *     .. External Functions ..
       LOGICAL            SISNAN, LSAME
       INTEGER            ICAMAX
-      REAL               SLAMCH
-      EXTERNAL           SISNAN, LSAME, ICAMAX, SLAMCH
+      REAL               SLAMCH, SCNRM2
+      EXTERNAL           SISNAN, LSAME, ICAMAX, SLAMCH, SCNRM2
 *     ..
 *     .. External Subroutines ..
       EXTERNAL           CSSCAL, CSWAP, XERBLA
@@ -325,15 +325,9 @@
       NOCONV = .FALSE.
 *
       DO 200 I = K, L
-         C = ZERO
-         R = ZERO
-*
-         DO 150 J = K, L
-            IF( J.EQ.I )
-     $         GO TO 150
-            C = C + CABS1( A( J, I ) )
-            R = R + CABS1( A( I, J ) )
-  150    CONTINUE
+*
+         C = SCNRM2( L-K+1, A( K, I ), 1 )
+         R = SCNRM2( L-K+1, A( I , K ), LDA )
          ICA = ICAMAX( L, A( 1, I ), 1 )
          CA = ABS( A( ICA, I ) )
          IRA = ICAMAX( N-K+1, A( I, K ), LDA )
diff --git a/SRC/dgebal.f b/SRC/dgebal.f
index 5d7ed035..fe4c156b 100644
--- a/SRC/dgebal.f
+++ b/SRC/dgebal.f
@@ -192,8 +192,8 @@
 *     .. External Functions ..
       LOGICAL            DISNAN, LSAME
       INTEGER            IDAMAX
-      DOUBLE PRECISION   DLAMCH
-      EXTERNAL           DISNAN, LSAME, IDAMAX, DLAMCH
+      DOUBLE PRECISION   DLAMCH, DNRM2
+      EXTERNAL           DISNAN, LSAME, IDAMAX, DLAMCH, DNRM2
 *     ..
 *     .. External Subroutines ..
       EXTERNAL           DSCAL, DSWAP, XERBLA
@@ -312,19 +312,14 @@
       SFMAX1 = ONE / SFMIN1
       SFMIN2 = SFMIN1*SCLFAC
       SFMAX2 = ONE / SFMIN2
+*
   140 CONTINUE
       NOCONV = .FALSE.
 *
       DO 200 I = K, L
-         C = ZERO
-         R = ZERO
-*
-         DO 150 J = K, L
-            IF( J.EQ.I )
-     $         GO TO 150
-            C = C + ABS( A( J, I ) )
-            R = R + ABS( A( I, J ) )
-  150    CONTINUE
+*
+         C = DNRM2( L-K+1, A( K, I ), 1 )
+         R = DNRM2( L-K+1, A( I, K ), LDA )
          ICA = IDAMAX( L, A( 1, I ), 1 )
          CA = ABS( A( ICA, I ) )
          IRA = IDAMAX( N-K+1, A( I, K ), LDA )
diff --git a/SRC/sgebal.f b/SRC/sgebal.f
index 853ff20e..de948488 100644
--- a/SRC/sgebal.f
+++ b/SRC/sgebal.f
@@ -192,8 +192,8 @@
 *     .. External Functions ..
       LOGICAL            SISNAN, LSAME
       INTEGER            ISAMAX
-      REAL               SLAMCH
-      EXTERNAL           SISNAN, LSAME, ISAMAX, SLAMCH
+      REAL               SLAMCH, SNRM2
+      EXTERNAL           SISNAN, LSAME, ISAMAX, SLAMCH, SNRM2
 *     ..
 *     .. External Subroutines ..
       EXTERNAL           SSCAL, SSWAP, XERBLA
@@ -316,15 +316,9 @@
       NOCONV = .FALSE.
 *
       DO 200 I = K, L
-         C = ZERO
-         R = ZERO
-*
-         DO 150 J = K, L
-            IF( J.EQ.I )
-     $         GO TO 150
-            C = C + ABS( A( J, I ) )
-            R = R + ABS( A( I, J ) )
-  150    CONTINUE
+*
+         C = SNRM2( L-K+1, A( K, I ), 1 )
+         R = SNRM2( L-K+1, A( I, K ), LDA )
          ICA = ISAMAX( L, A( 1, I ), 1 )
          CA = ABS( A( ICA, I ) )
          IRA = ISAMAX( N-K+1, A( I, K ), LDA )
diff --git a/SRC/zgebal.f b/SRC/zgebal.f
index 9c90f0b4..aa4f6d0e 100644
--- a/SRC/zgebal.f
+++ b/SRC/zgebal.f
@@ -194,8 +194,8 @@
 *     .. External Functions ..
       LOGICAL            DISNAN, LSAME
       INTEGER            IZAMAX
-      DOUBLE PRECISION   DLAMCH
-      EXTERNAL           DISNAN, LSAME, IZAMAX, DLAMCH
+      DOUBLE PRECISION   DLAMCH, DZNRM2
+      EXTERNAL           DISNAN, LSAME, IZAMAX, DLAMCH, DZNRM2
 *     ..
 *     .. External Subroutines ..
       EXTERNAL           XERBLA, ZDSCAL, ZSWAP
@@ -324,15 +324,9 @@
       NOCONV = .FALSE.
 *
       DO 200 I = K, L
-         C = ZERO
-         R = ZERO
-*
-         DO 150 J = K, L
-            IF( J.EQ.I )
-     $         GO TO 150
-            C = C + CABS1( A( J, I ) )
-            R = R + CABS1( A( I, J ) )
-  150    CONTINUE
+*
+         C = DZNRM2( L-K+1, A( K, I ), 1 )
+         R = DZNRM2( L-K+1, A( I, K ), LDA )
          ICA = IZAMAX( L, A( 1, I ), 1 )
          CA = ABS( A( ICA, I ) )
          IRA = IZAMAX( N-K+1, A( I, K ), LDA )