==============================================================================

Patch from Christof Voemel, ETH Zurich.
==============================================================================
I would like to make a small change  in dlarrd/slarrd/_larrd.f.base so that the
code deals better with certain matrices from Godunov.

Could you please replace in those files the lines, after label 40,
            GL = GL - FUDGE*SPDIAM*EPS*IN - FUDGE*PIVMIN
            GU = GU + FUDGE*SPDIAM*EPS*IN + FUDGE*PIVMIN

by
            GL = GL - FUDGE*TNORM*EPS*IN - FUDGE*PIVMIN
            GU = GU + FUDGE*TNORM*EPS*IN + FUDGE*PIVMIN

(Exchange SPDIAM to TNORM).
==============================================================================

diff --git a/SRC/dlarrd.f b/SRC/dlarrd.f
index 21effabf20eab7507f7b073577bd37e4307a4a06..d4b2e7ba58b094291fe593561fe8497922376b29 100644 (file)
--- a/SRC/dlarrd.f
+++ b/SRC/dlarrd.f
@@ -218,7 +218,7 @@
      $                   IM, IN, IOFF, IOUT, IRANGE, ITMAX, ITMP1,
      $                   ITMP2, IW, IWOFF, J, JBLK, JDISC, JE, JEE, NB,
      $                   NWL, NWU
-      DOUBLE PRECISION   ATOLI, EPS, GL, GU, RTOLI, SPDIAM, TMP1, TMP2,
+      DOUBLE PRECISION   ATOLI, EPS, GL, GU, RTOLI, TMP1, TMP2,
      $                   TNORM, UFLOW, WKILL, WLU, WUL
 
 *     ..
@@ -325,7 +325,8 @@
       TNORM = MAX( ABS( GL ), ABS( GU ) )
       GL = GL - FUDGE*TNORM*EPS*N - FUDGE*TWO*PIVMIN
       GU = GU + FUDGE*TNORM*EPS*N + FUDGE*TWO*PIVMIN
-      SPDIAM = GU - GL
+*     [JAN/28/2009] remove the line below since SPDIAM variable not use
+*     SPDIAM = GU - GL
 *     Input arguments for DLAEBZ:
 *     The relative tolerance.  An interval (a,b] lies within
 *     "relative tolerance" if  b-a < RELTOL*max(|a|,|b|),
@@ -490,9 +491,14 @@
                GL =  MIN( GL, GERS( 2*J - 1))
                GU = MAX( GU, GERS(2*J) )
    40       CONTINUE
-            SPDIAM = GU - GL
-            GL = GL - FUDGE*SPDIAM*EPS*IN - FUDGE*PIVMIN
-            GU = GU + FUDGE*SPDIAM*EPS*IN + FUDGE*PIVMIN
+*           [JAN/28/2009]
+*           change SPDIAM by TNORM in lines 2 and 3 thereafter
+*           line 1: remove computation of SPDIAM (not useful anymore)
+*           SPDIAM = GU - GL
+*           GL = GL - FUDGE*SPDIAM*EPS*IN - FUDGE*PIVMIN
+*           GU = GU + FUDGE*SPDIAM*EPS*IN + FUDGE*PIVMIN
+            GL = GL - FUDGE*TNORM*EPS*IN - FUDGE*PIVMIN
+            GU = GU + FUDGE*TNORM*EPS*IN + FUDGE*PIVMIN
 *
             IF( IRANGE.GT.1 ) THEN
                IF( GU.LT.WL ) THEN

diff --git a/SRC/slarrd.f b/SRC/slarrd.f
index 70c3ff03841a6ec8fa13dac732fed02855418d5a..1c63a2d994a849796b03c95bceb3e39581d2f862 100644 (file)
--- a/SRC/slarrd.f
+++ b/SRC/slarrd.f
@@ -218,7 +218,7 @@
      $                   IM, IN, IOFF, IOUT, IRANGE, ITMAX, ITMP1,
      $                   ITMP2, IW, IWOFF, J, JBLK, JDISC, JE, JEE, NB,
      $                   NWL, NWU
-      REAL               ATOLI, EPS, GL, GU, RTOLI, SPDIAM, TMP1, TMP2,
+      REAL               ATOLI, EPS, GL, GU, RTOLI, TMP1, TMP2,
      $                   TNORM, UFLOW, WKILL, WLU, WUL
 
 *     ..
@@ -325,7 +325,8 @@
       TNORM = MAX( ABS( GL ), ABS( GU ) )
       GL = GL - FUDGE*TNORM*EPS*N - FUDGE*TWO*PIVMIN
       GU = GU + FUDGE*TNORM*EPS*N + FUDGE*TWO*PIVMIN
-      SPDIAM = GU - GL
+*     [JAN/28/2009] remove the line below since SPDIAM variable not use
+*     SPDIAM = GU - GL
 *     Input arguments for SLAEBZ:
 *     The relative tolerance.  An interval (a,b] lies within
 *     "relative tolerance" if  b-a < RELTOL*max(|a|,|b|),
@@ -490,9 +491,14 @@
                GL =  MIN( GL, GERS( 2*J - 1))
                GU = MAX( GU, GERS(2*J) )
    40       CONTINUE
-            SPDIAM = GU - GL
-            GL = GL - FUDGE*SPDIAM*EPS*IN - FUDGE*PIVMIN
-            GU = GU + FUDGE*SPDIAM*EPS*IN + FUDGE*PIVMIN
+*           [JAN/28/2009]
+*           change SPDIAM by TNORM in lines 2 and 3 thereafter
+*           line 1: remove computation of SPDIAM (not useful anymore)
+*           SPDIAM = GU - GL
+*           GL = GL - FUDGE*SPDIAM*EPS*IN - FUDGE*PIVMIN
+*           GU = GU + FUDGE*SPDIAM*EPS*IN + FUDGE*PIVMIN
+            GL = GL - FUDGE*TNORM*EPS*IN - FUDGE*PIVMIN
+            GU = GU + FUDGE*TNORM*EPS*IN + FUDGE*PIVMIN
 *
             IF( IRANGE.GT.1 ) THEN
                IF( GU.LT.WL ) THEN