diff options
| author | Julien Langou <julien.langou@ucdenver.edu> | 2017-05-18 11:09:45 +0200 |
|---|---|---|
| committer | Julien Langou <julien.langou@ucdenver.edu> | 2017-05-18 11:09:45 +0200 |
| commit | a1e853bd005a2107a0a43f92782a158a9fd6e56c (patch) | |
| tree | 540f91580255a5911e61e5c7cbb2c7714690cb1c | |
| parent | 2f60afb4f081a02402e9b434cd6f151d847dd151 (diff) | |
| download | lapack-a1e853bd005a2107a0a43f92782a158a9fd6e56c.tar.gz lapack-a1e853bd005a2107a0a43f92782a158a9fd6e56c.tar.bz2 lapack-a1e853bd005a2107a0a43f92782a158a9fd6e56c.zip | |
Fix #147: xlapy2 not propagating nans
xLAPY2 now returns a NaN whenever input variables X or Y are NaNs.
The previous xLAPY2 was relying on FORTRAN INTRINSIC MAX and MIN to behave in a
certain way with NaNs (i.e. return a NaN whenever X or Y are NaN) and this
behavior is not observed on some (most?) compilers.
We handle the NaN behavior of xLAPY2 by checking for NaNs at the start of the
function.
Thanks to Andreas Noack for providing report and sample code to demonstrate the
problem.
| -rw-r--r-- | SRC/dlapy2.f | 28 | ||||
| -rw-r--r-- | SRC/slapy2.f | 31 |
2 files changed, 43 insertions, 16 deletions
diff --git a/SRC/dlapy2.f b/SRC/dlapy2.f index 3861b1d0..b40c46af 100644 --- a/SRC/dlapy2.f +++ b/SRC/dlapy2.f @@ -82,20 +82,32 @@ * .. * .. Local Scalars .. DOUBLE PRECISION W, XABS, YABS, Z + LOGICAL X_IS_NAN, Y_IS_NAN +* .. +* .. External Functions .. + LOGICAL DISNAN + EXTERNAL DISNAN * .. * .. Intrinsic Functions .. INTRINSIC ABS, MAX, MIN, SQRT * .. * .. Executable Statements .. * - XABS = ABS( X ) - YABS = ABS( Y ) - W = MAX( XABS, YABS ) - Z = MIN( XABS, YABS ) - IF( Z.EQ.ZERO ) THEN - DLAPY2 = W - ELSE - DLAPY2 = W*SQRT( ONE+( Z / W )**2 ) + X_IS_NAN = DISNAN( X ) + Y_IS_NAN = DISNAN( Y ) + IF ( X_IS_NAN ) DLAPY2 = X + IF ( Y_IS_NAN ) DLAPY2 = Y +* + IF ( .NOT.( X_IS_NAN.OR.Y_IS_NAN ) ) THEN + XABS = ABS( X ) + YABS = ABS( Y ) + W = MAX( XABS, YABS ) + Z = MIN( XABS, YABS ) + IF( Z.EQ.ZERO ) THEN + DLAPY2 = W + ELSE + DLAPY2 = W*SQRT( ONE+( Z / W )**2 ) + END IF END IF RETURN * diff --git a/SRC/slapy2.f b/SRC/slapy2.f index 13e21981..a8ee8b5b 100644 --- a/SRC/slapy2.f +++ b/SRC/slapy2.f @@ -82,20 +82,35 @@ * .. * .. Local Scalars .. REAL W, XABS, YABS, Z + LOGICAL X_IS_NAN, Y_IS_NAN +* .. +* .. External Functions .. + LOGICAL SISNAN + EXTERNAL SISNAN * .. * .. Intrinsic Functions .. INTRINSIC ABS, MAX, MIN, SQRT * .. * .. Executable Statements .. * - XABS = ABS( X ) - YABS = ABS( Y ) - W = MAX( XABS, YABS ) - Z = MIN( XABS, YABS ) - IF( Z.EQ.ZERO ) THEN - SLAPY2 = W - ELSE - SLAPY2 = W*SQRT( ONE+( Z / W )**2 ) +* .. +* .. Executable Statements .. +* + X_IS_NAN = SISNAN( X ) + Y_IS_NAN = SISNAN( Y ) + IF ( X_IS_NAN ) SLAPY2 = X + IF ( Y_IS_NAN ) SLAPY2 = Y +* + IF ( .NOT.( X_IS_NAN.OR.Y_IS_NAN ) ) THEN + XABS = ABS( X ) + YABS = ABS( Y ) + W = MAX( XABS, YABS ) + Z = MIN( XABS, YABS ) + IF( Z.EQ.ZERO ) THEN + SLAPY2 = W + ELSE + SLAPY2 = W*SQRT( ONE+( Z / W )**2 ) + END IF END IF RETURN * |