Applied patch provided by Sergey Kuznetsov @ Intel, fixes bug where TAU2 can be unini...

diff --git a/SRC/dlasd4.f b/SRC/dlasd4.f
index 877560c..838f716 100644 (file)
--- a/SRC/dlasd4.f
+++ b/SRC/dlasd4.f
@@ -223,6 +223,7 @@
 *
       EPS = DLAMCH( 'Epsilon' )
       RHOINV = ONE / RHO
+      TAU2= ZERO
 *
 *     The case I = N
 *
@@ -275,6 +276,7 @@
                ELSE
                   TAU2 = ( A+SQRT( A*A+FOUR*B*C ) ) / ( TWO*C )
                END IF
+               TAU = TAU2 / ( D( N )+SQRT( D( N )*D( N )+TAU2 ) )
             END IF
 *
 *           It can be proved that
@@ -293,6 +295,8 @@
             ELSE
                TAU2 = ( A+SQRT( A*A+FOUR*B*C ) ) / ( TWO*C )
             END IF
+            TAU = TAU2 / ( D( N )+SQRT( D( N )*D( N )+TAU2 ) )
+
 *
 *           It can be proved that
 *           D(N)^2 < D(N)^2+TAU2 < SIGMA(N)^2 < D(N)^2+RHO/2
@@ -301,7 +305,7 @@
 *
 *        The following TAU is to approximate SIGMA_n - D( N )
 *
-         TAU = TAU2 / ( D( N )+SQRT( D( N )*D( N )+TAU2 ) )
+*         TAU = TAU2 / ( D( N )+SQRT( D( N )*D( N )+TAU2 ) )
 *
          SIGMA = D( N ) + TAU
          DO 30 J = 1, N

diff --git a/SRC/slasd4.f b/SRC/slasd4.f
index 6bed144..bf0f2cf 100644 (file)
--- a/SRC/slasd4.f
+++ b/SRC/slasd4.f
@@ -223,6 +223,7 @@
 *
       EPS = SLAMCH( 'Epsilon' )
       RHOINV = ONE / RHO
+      TAU2= ZERO
 *
 *     The case I = N
 *
@@ -275,6 +276,7 @@
                ELSE
                   TAU2 = ( A+SQRT( A*A+FOUR*B*C ) ) / ( TWO*C )
                END IF
+               TAU = TAU2 / ( D( N )+SQRT( D( N )*D( N )+TAU2 ) )
             END IF
 *
 *           It can be proved that
@@ -293,6 +295,8 @@
             ELSE
                TAU2 = ( A+SQRT( A*A+FOUR*B*C ) ) / ( TWO*C )
             END IF
+            TAU = TAU2 / ( D( N )+SQRT( D( N )*D( N )+TAU2 ) )
+
 *
 *           It can be proved that
 *           D(N)^2 < D(N)^2+TAU2 < SIGMA(N)^2 < D(N)^2+RHO/2
@@ -301,7 +305,7 @@
 *
 *        The following TAU is to approximate SIGMA_n - D( N )
 *
-         TAU = TAU2 / ( D( N )+SQRT( D( N )*D( N )+TAU2 ) )
+*         TAU = TAU2 / ( D( N )+SQRT( D( N )*D( N )+TAU2 ) )
 *
          SIGMA = D( N ) + TAU
          DO 30 J = 1, N