summaryrefslogtreecommitdiff
path: root/src/math
diff options
context:
space:
mode:
authorSzabolcs Nagy <nsz@port70.net>2013-05-19 14:43:32 +0000
committerSzabolcs Nagy <nsz@port70.net>2013-05-19 14:43:32 +0000
commit1e5eb73545ca6cfe8b918798835aaf6e07af5beb (patch)
treed1af43fb6a86387586e5fbdcede7d2a3f56c78d4 /src/math
parentffd8ac2dd50f99c3c83d7d9d845df9874ec3e7d5 (diff)
downloadmusl-1e5eb73545ca6cfe8b918798835aaf6e07af5beb.tar.gz
math: add fma TODO comments about the underflow issue
The underflow exception is not raised correctly in some cornercases (see previous fma commit), added comments with examples for fmaf, fmal and non-x86 fma. In fmaf store the result before returning so it has the correct precision when FLT_EVAL_METHOD!=0
Diffstat (limited to 'src/math')
-rw-r--r--src/math/fma.c2
-rw-r--r--src/math/fmaf.c12
-rw-r--r--src/math/fmal.c2
3 files changed, 14 insertions, 2 deletions
diff --git a/src/math/fma.c b/src/math/fma.c
index 89def795..850a4a6c 100644
--- a/src/math/fma.c
+++ b/src/math/fma.c
@@ -441,6 +441,8 @@ double fma(double x, double y, double z)
/*
* There is no need to worry about double rounding in directed
* rounding modes.
+ * TODO: underflow is not raised properly, example in downward rounding:
+ * fma(0x1.000000001p-1000, 0x1.000000001p-30, -0x1p-1066)
*/
fesetround(oround);
adj = r.lo + xy.lo;
diff --git a/src/math/fmaf.c b/src/math/fmaf.c
index a1c7f4f8..745ee393 100644
--- a/src/math/fmaf.c
+++ b/src/math/fmaf.c
@@ -49,7 +49,14 @@ float fmaf(float x, float y, float z)
(hr & 0x7ff00000) == 0x7ff00000 || /* NaN */
result - xy == z || /* exact */
fegetround() != FE_TONEAREST) /* not round-to-nearest */
- return (result);
+ {
+ /*
+ TODO: underflow is not raised correctly, example in
+ downward rouding: fmaf(0x1p-120f, 0x1p-120f, 0x1p-149f)
+ */
+ z = result;
+ return z;
+ }
/*
* If result is inexact, and exactly halfway between two float values,
@@ -63,5 +70,6 @@ float fmaf(float x, float y, float z)
fesetround(FE_TONEAREST);
if (result == adjusted_result)
SET_LOW_WORD(adjusted_result, lr + 1);
- return (adjusted_result);
+ z = adjusted_result;
+ return z;
}
diff --git a/src/math/fmal.c b/src/math/fmal.c
index ccbe434d..87e30fcf 100644
--- a/src/math/fmal.c
+++ b/src/math/fmal.c
@@ -262,6 +262,8 @@ long double fmal(long double x, long double y, long double z)
/*
* There is no need to worry about double rounding in directed
* rounding modes.
+ * TODO: underflow is not raised correctly, example in downward rounding:
+ * fmal(0x1.0000000001p-16000L, 0x1.0000000001p-400L, -0x1p-16440L)
*/
fesetround(oround);
adj = r.lo + xy.lo;