From f4b5da4ee2f9d8505129dde7a2b1c01c3025e919 Mon Sep 17 00:00:00 2001 From: Simon Maertens Date: Thu, 10 Sep 2026 09:14:41 +0100 Subject: [PATCH 1/2] xLAEIN: apply the growth test to the growth, not to |x| alone The residual of the vector a solve returns is SCALE*|v|/|x|, so the acceptance test belongs on that ratio. Testing |x| against GROWTO alone made it depend on the starting vector: GROWTO suits the restart vectors, whose 1-norm is of order SQRT(N)*EPS3, and was SQRT(N) too lenient for the initial vector, whose 1-norm is N*EPS3 - and that is where almost every vector is accepted. One NEP matrix reached 27.2 against a test threshold of 20; the worst ratio is now 4.2. Co-Authored-By: Claude Opus 5 (1M context) --- SRC/claein.f | 21 +++++++++++++-------- SRC/dlaein.f | 35 ++++++++++++++++++++--------------- SRC/slaein.f | 35 ++++++++++++++++++++--------------- SRC/zlaein.f | 21 +++++++++++++-------- 4 files changed, 66 insertions(+), 46 deletions(-) diff --git a/SRC/claein.f b/SRC/claein.f index efc659a95..01651cb59 100644 --- a/SRC/claein.f +++ b/SRC/claein.f @@ -173,7 +173,8 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * .. Local Scalars .. CHARACTER NORMIN, TRANS INTEGER I, IERR, ITS, J - REAL GROWTO, NRMSML, ROOTN, RTEMP, SCALE, VNORM + REAL GROWTO, NRMSML, ROOTN, RTEMP, SCALE, VI_NORM, + $ VIP1_NORM COMPLEX CDUM, EI, EJ, TEMP, X * .. * .. External Functions .. @@ -198,11 +199,14 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * INFO = 0 * -* GROWTO is the threshold used in the acceptance test for an -* eigenvector. +* The residual of the vector x that a solve returns is SCALE times +* the norm of the starting vector, over the norm of x, so GROWTO is +* the growth VIP1_NORM/(SCALE*VI_NORM) that the acceptance test +* below requires, where VI_NORM and VIP1_NORM are the norms of the +* vectors the current iteration started from and produced. * ROOTN = SQRT( REAL( N ) ) - GROWTO = TENTH / ROOTN + GROWTO = TENTH / ( REAL( N )*EPS3 ) NRMSML = MAX( ONE, EPS3*ROOTN )*SMLNUM * * Form B = H - W*I (except that the subdiagonal elements are not @@ -226,8 +230,8 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * * Scale supplied initial vector. * - VNORM = SCNRM2( N, V, 1 ) - CALL CSSCAL( N, ( EPS3*ROOTN ) / MAX( VNORM, NRMSML ), V, + VI_NORM = SCNRM2( N, V, 1 ) + CALL CSSCAL( N, ( EPS3*ROOTN ) / MAX( VI_NORM, NRMSML ), V, $ 1 ) END IF * @@ -309,6 +313,7 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * NORMIN = 'N' DO 110 ITS = 1, N + VI_NORM = SCASUM( N, V, 1 ) * * Solve U*x = scale*v for a right eigenvector * or U**H *x = scale*v for a left eigenvector, @@ -321,8 +326,8 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * * Test for sufficient growth in the norm of v. * - VNORM = SCASUM( N, V, 1 ) - IF( VNORM.GE.GROWTO*SCALE ) + VIP1_NORM = SCASUM( N, V, 1 ) + IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) $ GO TO 120 * * Choose new orthogonal starting vector and try again. diff --git a/SRC/dlaein.f b/SRC/dlaein.f index 17bb819aa..c092def4f 100644 --- a/SRC/dlaein.f +++ b/SRC/dlaein.f @@ -194,8 +194,8 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, CHARACTER NORMIN, TRANS INTEGER I, I1, I2, I3, IERR, ITS, J DOUBLE PRECISION ABSBII, ABSBJJ, EI, EJ, GROWTO, NORM, NRMSML, - $ REC, ROOTN, SCALE, TEMP, VCRIT, VMAX, VNORM, W, - $ W1, X, XI, XR, Y + $ REC, ROOTN, SCALE, TEMP, VCRIT, VI_NORM, + $ VIP1_NORM, VMAX, W, W1, X, XI, XR, Y * .. * .. External Functions .. INTEGER IDAMAX @@ -212,11 +212,14 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * INFO = 0 * -* GROWTO is the threshold used in the acceptance test for an -* eigenvector. +* The residual of the vector x that a solve returns is SCALE times +* the norm of the starting vector, over the norm of x, so GROWTO is +* the growth VIP1_NORM/(SCALE*VI_NORM) that the acceptance test +* below requires, where VI_NORM and VIP1_NORM are the norms of the +* vectors the current iteration started from and produced. * ROOTN = SQRT( DBLE( N ) ) - GROWTO = TENTH / ROOTN + GROWTO = TENTH / ( DBLE( N )*EPS3 ) NRMSML = MAX( ONE, EPS3*ROOTN )*SMLNUM * * Form B = H - (WR,WI)*I (except that the subdiagonal elements and @@ -244,8 +247,8 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Scale supplied initial vector. * - VNORM = DNRM2( N, VR, 1 ) - CALL DSCAL( N, ( EPS3*ROOTN ) / MAX( VNORM, NRMSML ), VR, + VI_NORM = DNRM2( N, VR, 1 ) + CALL DSCAL( N, ( EPS3*ROOTN ) / MAX( VI_NORM, NRMSML ), VR, $ 1 ) END IF * @@ -327,6 +330,7 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * NORMIN = 'N' DO 110 ITS = 1, N + VI_NORM = DASUM( N, VR, 1 ) * * Solve U*x = scale*v for a right eigenvector * or U**T*x = scale*v for a left eigenvector, @@ -339,8 +343,8 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Test for sufficient growth in the norm of v. * - VNORM = DASUM( N, VR, 1 ) - IF( VNORM.GE.GROWTO*SCALE ) + VIP1_NORM = DASUM( N, VR, 1 ) + IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) $ GO TO 120 * * Choose new orthogonal starting vector and try again. @@ -521,6 +525,7 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, END IF * DO 270 ITS = 1, N + VI_NORM = DASUM( N, VR, 1 ) + DASUM( N, VI, 1 ) SCALE = ONE VMAX = ONE VCRIT = BIGNUM @@ -591,8 +596,8 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Test for sufficient growth in the norm of (VR,VI). * - VNORM = DASUM( N, VR, 1 ) + DASUM( N, VI, 1 ) - IF( VNORM.GE.GROWTO*SCALE ) + VIP1_NORM = DASUM( N, VR, 1 ) + DASUM( N, VI, 1 ) + IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) $ GO TO 280 * * Choose a new orthogonal starting vector and try again. @@ -616,12 +621,12 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Normalize eigenvector. * - VNORM = ZERO + VIP1_NORM = ZERO DO 290 I = 1, N - VNORM = MAX( VNORM, ABS( VR( I ) )+ABS( VI( I ) ) ) + VIP1_NORM = MAX( VIP1_NORM, ABS( VR( I ) )+ABS( VI( I ) ) ) 290 CONTINUE - CALL DSCAL( N, ONE / VNORM, VR, 1 ) - CALL DSCAL( N, ONE / VNORM, VI, 1 ) + CALL DSCAL( N, ONE / VIP1_NORM, VR, 1 ) + CALL DSCAL( N, ONE / VIP1_NORM, VI, 1 ) * END IF * diff --git a/SRC/slaein.f b/SRC/slaein.f index e6e2065f5..93a168960 100644 --- a/SRC/slaein.f +++ b/SRC/slaein.f @@ -194,8 +194,8 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, CHARACTER NORMIN, TRANS INTEGER I, I1, I2, I3, IERR, ITS, J REAL ABSBII, ABSBJJ, EI, EJ, GROWTO, NORM, NRMSML, - $ REC, ROOTN, SCALE, TEMP, VCRIT, VMAX, VNORM, W, - $ W1, X, XI, XR, Y + $ REC, ROOTN, SCALE, TEMP, VCRIT, VI_NORM, + $ VIP1_NORM, VMAX, W, W1, X, XI, XR, Y * .. * .. External Functions .. INTEGER ISAMAX @@ -212,11 +212,14 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * INFO = 0 * -* GROWTO is the threshold used in the acceptance test for an -* eigenvector. +* The residual of the vector x that a solve returns is SCALE times +* the norm of the starting vector, over the norm of x, so GROWTO is +* the growth VIP1_NORM/(SCALE*VI_NORM) that the acceptance test +* below requires, where VI_NORM and VIP1_NORM are the norms of the +* vectors the current iteration started from and produced. * ROOTN = SQRT( REAL( N ) ) - GROWTO = TENTH / ROOTN + GROWTO = TENTH / ( REAL( N )*EPS3 ) NRMSML = MAX( ONE, EPS3*ROOTN )*SMLNUM * * Form B = H - (WR,WI)*I (except that the subdiagonal elements and @@ -244,8 +247,8 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Scale supplied initial vector. * - VNORM = SNRM2( N, VR, 1 ) - CALL SSCAL( N, ( EPS3*ROOTN ) / MAX( VNORM, NRMSML ), VR, + VI_NORM = SNRM2( N, VR, 1 ) + CALL SSCAL( N, ( EPS3*ROOTN ) / MAX( VI_NORM, NRMSML ), VR, $ 1 ) END IF * @@ -327,6 +330,7 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * NORMIN = 'N' DO 110 ITS = 1, N + VI_NORM = SASUM( N, VR, 1 ) * * Solve U*x = scale*v for a right eigenvector * or U**T*x = scale*v for a left eigenvector, @@ -339,8 +343,8 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Test for sufficient growth in the norm of v. * - VNORM = SASUM( N, VR, 1 ) - IF( VNORM.GE.GROWTO*SCALE ) + VIP1_NORM = SASUM( N, VR, 1 ) + IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) $ GO TO 120 * * Choose new orthogonal starting vector and try again. @@ -521,6 +525,7 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, END IF * DO 270 ITS = 1, N + VI_NORM = SASUM( N, VR, 1 ) + SASUM( N, VI, 1 ) SCALE = ONE VMAX = ONE VCRIT = BIGNUM @@ -591,8 +596,8 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Test for sufficient growth in the norm of (VR,VI). * - VNORM = SASUM( N, VR, 1 ) + SASUM( N, VI, 1 ) - IF( VNORM.GE.GROWTO*SCALE ) + VIP1_NORM = SASUM( N, VR, 1 ) + SASUM( N, VI, 1 ) + IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) $ GO TO 280 * * Choose a new orthogonal starting vector and try again. @@ -616,12 +621,12 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Normalize eigenvector. * - VNORM = ZERO + VIP1_NORM = ZERO DO 290 I = 1, N - VNORM = MAX( VNORM, ABS( VR( I ) )+ABS( VI( I ) ) ) + VIP1_NORM = MAX( VIP1_NORM, ABS( VR( I ) )+ABS( VI( I ) ) ) 290 CONTINUE - CALL SSCAL( N, ONE / VNORM, VR, 1 ) - CALL SSCAL( N, ONE / VNORM, VI, 1 ) + CALL SSCAL( N, ONE / VIP1_NORM, VR, 1 ) + CALL SSCAL( N, ONE / VIP1_NORM, VI, 1 ) * END IF * diff --git a/SRC/zlaein.f b/SRC/zlaein.f index 78b72b881..744007b78 100644 --- a/SRC/zlaein.f +++ b/SRC/zlaein.f @@ -173,7 +173,8 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * .. Local Scalars .. CHARACTER NORMIN, TRANS INTEGER I, IERR, ITS, J - DOUBLE PRECISION GROWTO, NRMSML, ROOTN, RTEMP, SCALE, VNORM + DOUBLE PRECISION GROWTO, NRMSML, ROOTN, RTEMP, SCALE, VI_NORM, + $ VIP1_NORM COMPLEX*16 CDUM, EI, EJ, TEMP, X * .. * .. External Functions .. @@ -198,11 +199,14 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * INFO = 0 * -* GROWTO is the threshold used in the acceptance test for an -* eigenvector. +* The residual of the vector x that a solve returns is SCALE times +* the norm of the starting vector, over the norm of x, so GROWTO is +* the growth VIP1_NORM/(SCALE*VI_NORM) that the acceptance test +* below requires, where VI_NORM and VIP1_NORM are the norms of the +* vectors the current iteration started from and produced. * ROOTN = SQRT( DBLE( N ) ) - GROWTO = TENTH / ROOTN + GROWTO = TENTH / ( DBLE( N )*EPS3 ) NRMSML = MAX( ONE, EPS3*ROOTN )*SMLNUM * * Form B = H - W*I (except that the subdiagonal elements are not @@ -226,8 +230,8 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * * Scale supplied initial vector. * - VNORM = DZNRM2( N, V, 1 ) - CALL ZDSCAL( N, ( EPS3*ROOTN ) / MAX( VNORM, NRMSML ), V, + VI_NORM = DZNRM2( N, V, 1 ) + CALL ZDSCAL( N, ( EPS3*ROOTN ) / MAX( VI_NORM, NRMSML ), V, $ 1 ) END IF * @@ -309,6 +313,7 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * NORMIN = 'N' DO 110 ITS = 1, N + VI_NORM = DZASUM( N, V, 1 ) * * Solve U*x = scale*v for a right eigenvector * or U**H *x = scale*v for a left eigenvector, @@ -321,8 +326,8 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * * Test for sufficient growth in the norm of v. * - VNORM = DZASUM( N, V, 1 ) - IF( VNORM.GE.GROWTO*SCALE ) + VIP1_NORM = DZASUM( N, V, 1 ) + IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) $ GO TO 120 * * Choose new orthogonal starting vector and try again. From dc0888a3f7f2593e491fdae62748c69345d96bb5 Mon Sep 17 00:00:00 2001 From: Simon Maertens Date: Thu, 10 Sep 2026 17:06:31 +0100 Subject: [PATCH 2/2] xLAEIN: call the two norms V0NORM and V1NORM Each starting vector gets exactly one solve, so there is no iterate beyond the first; 0 and 1 say that where i and i+1 do not. Co-Authored-By: Claude Opus 5 (1M context) --- SRC/claein.f | 23 +++++++++++------------ SRC/dlaein.f | 37 ++++++++++++++++++------------------- SRC/slaein.f | 37 ++++++++++++++++++------------------- SRC/zlaein.f | 23 +++++++++++------------ 4 files changed, 58 insertions(+), 62 deletions(-) diff --git a/SRC/claein.f b/SRC/claein.f index 01651cb59..bbd45755c 100644 --- a/SRC/claein.f +++ b/SRC/claein.f @@ -173,8 +173,8 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * .. Local Scalars .. CHARACTER NORMIN, TRANS INTEGER I, IERR, ITS, J - REAL GROWTO, NRMSML, ROOTN, RTEMP, SCALE, VI_NORM, - $ VIP1_NORM + REAL GROWTO, NRMSML, ROOTN, RTEMP, SCALE, V0NORM, + $ V1NORM COMPLEX CDUM, EI, EJ, TEMP, X * .. * .. External Functions .. @@ -199,11 +199,10 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * INFO = 0 * -* The residual of the vector x that a solve returns is SCALE times -* the norm of the starting vector, over the norm of x, so GROWTO is -* the growth VIP1_NORM/(SCALE*VI_NORM) that the acceptance test -* below requires, where VI_NORM and VIP1_NORM are the norms of the -* vectors the current iteration started from and produced. +* Each starting vector gets one solve, from v0 to v1. The residual +* of v1 is SCALE times the norm of v0, over the norm of v1, so +* GROWTO is the growth V1NORM/(SCALE*V0NORM) that the acceptance +* test below requires. * ROOTN = SQRT( REAL( N ) ) GROWTO = TENTH / ( REAL( N )*EPS3 ) @@ -230,8 +229,8 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * * Scale supplied initial vector. * - VI_NORM = SCNRM2( N, V, 1 ) - CALL CSSCAL( N, ( EPS3*ROOTN ) / MAX( VI_NORM, NRMSML ), V, + V0NORM = SCNRM2( N, V, 1 ) + CALL CSSCAL( N, ( EPS3*ROOTN ) / MAX( V0NORM, NRMSML ), V, $ 1 ) END IF * @@ -313,7 +312,7 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * NORMIN = 'N' DO 110 ITS = 1, N - VI_NORM = SCASUM( N, V, 1 ) + V0NORM = SCASUM( N, V, 1 ) * * Solve U*x = scale*v for a right eigenvector * or U**H *x = scale*v for a left eigenvector, @@ -326,8 +325,8 @@ SUBROUTINE CLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * * Test for sufficient growth in the norm of v. * - VIP1_NORM = SCASUM( N, V, 1 ) - IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) + V1NORM = SCASUM( N, V, 1 ) + IF( V1NORM.GE.GROWTO*SCALE*V0NORM ) $ GO TO 120 * * Choose new orthogonal starting vector and try again. diff --git a/SRC/dlaein.f b/SRC/dlaein.f index c092def4f..fbcab9b38 100644 --- a/SRC/dlaein.f +++ b/SRC/dlaein.f @@ -194,8 +194,8 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, CHARACTER NORMIN, TRANS INTEGER I, I1, I2, I3, IERR, ITS, J DOUBLE PRECISION ABSBII, ABSBJJ, EI, EJ, GROWTO, NORM, NRMSML, - $ REC, ROOTN, SCALE, TEMP, VCRIT, VI_NORM, - $ VIP1_NORM, VMAX, W, W1, X, XI, XR, Y + $ REC, ROOTN, SCALE, TEMP, V0NORM, V1NORM, + $ VCRIT, VMAX, W, W1, X, XI, XR, Y * .. * .. External Functions .. INTEGER IDAMAX @@ -212,11 +212,10 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * INFO = 0 * -* The residual of the vector x that a solve returns is SCALE times -* the norm of the starting vector, over the norm of x, so GROWTO is -* the growth VIP1_NORM/(SCALE*VI_NORM) that the acceptance test -* below requires, where VI_NORM and VIP1_NORM are the norms of the -* vectors the current iteration started from and produced. +* Each starting vector gets one solve, from v0 to v1. The residual +* of v1 is SCALE times the norm of v0, over the norm of v1, so +* GROWTO is the growth V1NORM/(SCALE*V0NORM) that the acceptance +* test below requires. * ROOTN = SQRT( DBLE( N ) ) GROWTO = TENTH / ( DBLE( N )*EPS3 ) @@ -247,8 +246,8 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Scale supplied initial vector. * - VI_NORM = DNRM2( N, VR, 1 ) - CALL DSCAL( N, ( EPS3*ROOTN ) / MAX( VI_NORM, NRMSML ), VR, + V0NORM = DNRM2( N, VR, 1 ) + CALL DSCAL( N, ( EPS3*ROOTN ) / MAX( V0NORM, NRMSML ), VR, $ 1 ) END IF * @@ -330,7 +329,7 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * NORMIN = 'N' DO 110 ITS = 1, N - VI_NORM = DASUM( N, VR, 1 ) + V0NORM = DASUM( N, VR, 1 ) * * Solve U*x = scale*v for a right eigenvector * or U**T*x = scale*v for a left eigenvector, @@ -343,8 +342,8 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Test for sufficient growth in the norm of v. * - VIP1_NORM = DASUM( N, VR, 1 ) - IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) + V1NORM = DASUM( N, VR, 1 ) + IF( V1NORM.GE.GROWTO*SCALE*V0NORM ) $ GO TO 120 * * Choose new orthogonal starting vector and try again. @@ -525,7 +524,7 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, END IF * DO 270 ITS = 1, N - VI_NORM = DASUM( N, VR, 1 ) + DASUM( N, VI, 1 ) + V0NORM = DASUM( N, VR, 1 ) + DASUM( N, VI, 1 ) SCALE = ONE VMAX = ONE VCRIT = BIGNUM @@ -596,8 +595,8 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Test for sufficient growth in the norm of (VR,VI). * - VIP1_NORM = DASUM( N, VR, 1 ) + DASUM( N, VI, 1 ) - IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) + V1NORM = DASUM( N, VR, 1 ) + DASUM( N, VI, 1 ) + IF( V1NORM.GE.GROWTO*SCALE*V0NORM ) $ GO TO 280 * * Choose a new orthogonal starting vector and try again. @@ -621,12 +620,12 @@ SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Normalize eigenvector. * - VIP1_NORM = ZERO + V1NORM = ZERO DO 290 I = 1, N - VIP1_NORM = MAX( VIP1_NORM, ABS( VR( I ) )+ABS( VI( I ) ) ) + V1NORM = MAX( V1NORM, ABS( VR( I ) )+ABS( VI( I ) ) ) 290 CONTINUE - CALL DSCAL( N, ONE / VIP1_NORM, VR, 1 ) - CALL DSCAL( N, ONE / VIP1_NORM, VI, 1 ) + CALL DSCAL( N, ONE / V1NORM, VR, 1 ) + CALL DSCAL( N, ONE / V1NORM, VI, 1 ) * END IF * diff --git a/SRC/slaein.f b/SRC/slaein.f index 93a168960..8246aeaf3 100644 --- a/SRC/slaein.f +++ b/SRC/slaein.f @@ -194,8 +194,8 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, CHARACTER NORMIN, TRANS INTEGER I, I1, I2, I3, IERR, ITS, J REAL ABSBII, ABSBJJ, EI, EJ, GROWTO, NORM, NRMSML, - $ REC, ROOTN, SCALE, TEMP, VCRIT, VI_NORM, - $ VIP1_NORM, VMAX, W, W1, X, XI, XR, Y + $ REC, ROOTN, SCALE, TEMP, V0NORM, V1NORM, + $ VCRIT, VMAX, W, W1, X, XI, XR, Y * .. * .. External Functions .. INTEGER ISAMAX @@ -212,11 +212,10 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * INFO = 0 * -* The residual of the vector x that a solve returns is SCALE times -* the norm of the starting vector, over the norm of x, so GROWTO is -* the growth VIP1_NORM/(SCALE*VI_NORM) that the acceptance test -* below requires, where VI_NORM and VIP1_NORM are the norms of the -* vectors the current iteration started from and produced. +* Each starting vector gets one solve, from v0 to v1. The residual +* of v1 is SCALE times the norm of v0, over the norm of v1, so +* GROWTO is the growth V1NORM/(SCALE*V0NORM) that the acceptance +* test below requires. * ROOTN = SQRT( REAL( N ) ) GROWTO = TENTH / ( REAL( N )*EPS3 ) @@ -247,8 +246,8 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Scale supplied initial vector. * - VI_NORM = SNRM2( N, VR, 1 ) - CALL SSCAL( N, ( EPS3*ROOTN ) / MAX( VI_NORM, NRMSML ), VR, + V0NORM = SNRM2( N, VR, 1 ) + CALL SSCAL( N, ( EPS3*ROOTN ) / MAX( V0NORM, NRMSML ), VR, $ 1 ) END IF * @@ -330,7 +329,7 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * NORMIN = 'N' DO 110 ITS = 1, N - VI_NORM = SASUM( N, VR, 1 ) + V0NORM = SASUM( N, VR, 1 ) * * Solve U*x = scale*v for a right eigenvector * or U**T*x = scale*v for a left eigenvector, @@ -343,8 +342,8 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Test for sufficient growth in the norm of v. * - VIP1_NORM = SASUM( N, VR, 1 ) - IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) + V1NORM = SASUM( N, VR, 1 ) + IF( V1NORM.GE.GROWTO*SCALE*V0NORM ) $ GO TO 120 * * Choose new orthogonal starting vector and try again. @@ -525,7 +524,7 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, END IF * DO 270 ITS = 1, N - VI_NORM = SASUM( N, VR, 1 ) + SASUM( N, VI, 1 ) + V0NORM = SASUM( N, VR, 1 ) + SASUM( N, VI, 1 ) SCALE = ONE VMAX = ONE VCRIT = BIGNUM @@ -596,8 +595,8 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Test for sufficient growth in the norm of (VR,VI). * - VIP1_NORM = SASUM( N, VR, 1 ) + SASUM( N, VI, 1 ) - IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) + V1NORM = SASUM( N, VR, 1 ) + SASUM( N, VI, 1 ) + IF( V1NORM.GE.GROWTO*SCALE*V0NORM ) $ GO TO 280 * * Choose a new orthogonal starting vector and try again. @@ -621,12 +620,12 @@ SUBROUTINE SLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, * * Normalize eigenvector. * - VIP1_NORM = ZERO + V1NORM = ZERO DO 290 I = 1, N - VIP1_NORM = MAX( VIP1_NORM, ABS( VR( I ) )+ABS( VI( I ) ) ) + V1NORM = MAX( V1NORM, ABS( VR( I ) )+ABS( VI( I ) ) ) 290 CONTINUE - CALL SSCAL( N, ONE / VIP1_NORM, VR, 1 ) - CALL SSCAL( N, ONE / VIP1_NORM, VI, 1 ) + CALL SSCAL( N, ONE / V1NORM, VR, 1 ) + CALL SSCAL( N, ONE / V1NORM, VI, 1 ) * END IF * diff --git a/SRC/zlaein.f b/SRC/zlaein.f index 744007b78..0786fda9f 100644 --- a/SRC/zlaein.f +++ b/SRC/zlaein.f @@ -173,8 +173,8 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * .. Local Scalars .. CHARACTER NORMIN, TRANS INTEGER I, IERR, ITS, J - DOUBLE PRECISION GROWTO, NRMSML, ROOTN, RTEMP, SCALE, VI_NORM, - $ VIP1_NORM + DOUBLE PRECISION GROWTO, NRMSML, ROOTN, RTEMP, SCALE, V0NORM, + $ V1NORM COMPLEX*16 CDUM, EI, EJ, TEMP, X * .. * .. External Functions .. @@ -199,11 +199,10 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * INFO = 0 * -* The residual of the vector x that a solve returns is SCALE times -* the norm of the starting vector, over the norm of x, so GROWTO is -* the growth VIP1_NORM/(SCALE*VI_NORM) that the acceptance test -* below requires, where VI_NORM and VIP1_NORM are the norms of the -* vectors the current iteration started from and produced. +* Each starting vector gets one solve, from v0 to v1. The residual +* of v1 is SCALE times the norm of v0, over the norm of v1, so +* GROWTO is the growth V1NORM/(SCALE*V0NORM) that the acceptance +* test below requires. * ROOTN = SQRT( DBLE( N ) ) GROWTO = TENTH / ( DBLE( N )*EPS3 ) @@ -230,8 +229,8 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * * Scale supplied initial vector. * - VI_NORM = DZNRM2( N, V, 1 ) - CALL ZDSCAL( N, ( EPS3*ROOTN ) / MAX( VI_NORM, NRMSML ), V, + V0NORM = DZNRM2( N, V, 1 ) + CALL ZDSCAL( N, ( EPS3*ROOTN ) / MAX( V0NORM, NRMSML ), V, $ 1 ) END IF * @@ -313,7 +312,7 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * NORMIN = 'N' DO 110 ITS = 1, N - VI_NORM = DZASUM( N, V, 1 ) + V0NORM = DZASUM( N, V, 1 ) * * Solve U*x = scale*v for a right eigenvector * or U**H *x = scale*v for a left eigenvector, @@ -326,8 +325,8 @@ SUBROUTINE ZLAEIN( RIGHTV, NOINIT, N, H, LDH, W, V, B, LDB, * * Test for sufficient growth in the norm of v. * - VIP1_NORM = DZASUM( N, V, 1 ) - IF( VIP1_NORM.GE.GROWTO*SCALE*VI_NORM ) + V1NORM = DZASUM( N, V, 1 ) + IF( V1NORM.GE.GROWTO*SCALE*V0NORM ) $ GO TO 120 * * Choose new orthogonal starting vector and try again.