Patch LAPACK XLASD4.f as discussed in JuliaLang/julia#2340

This commit is contained in:
Viral B. Shah 2013-08-21 19:14:07 +05:30
parent c92ae012a6
commit eae6920f2d
2 changed files with 10 additions and 2 deletions

View File

@ -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

View File

@ -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