Thread: exp() versus the POSIX standard

The POSIX standard says this about the exp(3) function:

    If x is -Inf, +0 shall be returned.

At least on my Linux box, our version does no such thing:

regression=# select exp('-inf'::float8);
ERROR:  value out of range: underflow

Does anyone disagree that that's a bug?  Should we back-patch
a fix, or just change it in HEAD?  Given the lack of user
complaints, I lean a bit towards the latter, but am not sure.

            regards, tom lane

I wrote:
> The POSIX standard says this about the exp(3) function:
>     If x is -Inf, +0 shall be returned.
> At least on my Linux box, our version does no such thing:
> regression=# select exp('-inf'::float8);
> ERROR:  value out of range: underflow

Now that I look, power() has similar issues:

regression=# select power('1.1'::float8, '-inf');
ERROR:  value out of range: underflow
regression=# select power('0.1'::float8, 'inf');
ERROR:  value out of range: underflow
regression=# select power('-inf'::float8, '-3');
ERROR:  value out of range: underflow
regression=# select power('-inf'::float8, '-4');
ERROR:  value out of range: underflow

contradicting POSIX which says

For |x| > 1, if y is -Inf, +0 shall be returned.

For |x| < 1, if y is +Inf, +0 shall be returned.

For y an odd integer < 0, if x is -Inf, -0 shall be returned.

For y < 0 and not an odd integer, if x is -Inf, +0 shall be returned.

            regards, tom lane

=?UTF-8?Q?Darafei_=22Kom=D1=8Fpa=22_Praliaskouski?= <me@komzpa.net> writes:
> I've had the same issue with multiplying two tiny numbers. Select
> 2e-300::float * 2e-300::float gives an underflow, and it is not a wanted
> thing. This looks like handmade implementation of IEEE754's underflow
> exception that should be an optional return flag in addition to well
> defined number, but became a stop-the-world exception instead.

Solving that problem is very far outside the scope of what I'm interested
in here.  I think that we'd probably regret it if we try to support IEEE
subnormals, for example --- I know that all modern hardware is probably
good with those, but I'd bet against different platforms' libc functions
all behaving the same.  I don't see a sane way to offer user control over
whether we throw underflow errors or not, either.  (Do you really want "+"
to stop being immutable?)  The darker corners of IEEE754, like inexactness
exceptions, are even less likely to be implemented consistently
everywhere.

            regards, tom lane

> Does anyone disagree that that's a bug?  Should we back-patch
> a fix, or just change it in HEAD?  Given the lack of user
> complaints, I lean a bit towards the latter, but am not sure.

The other functions and operators pay attention to not give an error
when the input is Inf or 0.   exp() and power() are at least
inconsistent by doing so.  I don't think this behavior is useful.
Although it'd still be less risky to fix it in HEAD only.

> No, no kind of GUC switch is needed. Just drop underflow/overflow checks. You'll get 0 or Infinity in expected
places,and infinities are okay since 2007.
 

This is out of scope of this thread.  I am not sure opening it to
discussion on another thread would yield any result.  Experienced
developers like Tom appear to be in agreement of us needing to protect
users from oddities of floating point numbers.  (I am not.)

Emre Hasegeli <emre@hasegeli.com> writes:
>> No, no kind of GUC switch is needed. Just drop underflow/overflow checks. You'll get 0 or Infinity in expected
places,and infinities are okay since 2007. 

> This is out of scope of this thread.

Yeah, that.  At the moment I'm just interested in making the float and
numeric functions give equivalent results for infinite inputs.  If you
want to make a more general proposal about removing error checks, that
seems like a separate topic.

> I am not sure opening it to
> discussion on another thread would yield any result.  Experienced
> developers like Tom appear to be in agreement of us needing to protect
> users from oddities of floating point numbers.  (I am not.)

I think there's a pretty fundamental distinction between this behavior:

regression=# select exp('-inf'::float8);
 exp
-----
   0
(1 row)

and this one:

regression=# select exp('-1000'::float8);
ERROR:  value out of range: underflow

In the first case, zero is the correct answer to any precision you care
to name.  In the second case, zero is *not* the correct answer; we simply
cannot represent the correct answer (somewhere around 1e-434) as a float8.
Returning zero would represent 100% loss of accuracy.  Now, there may well
be applications where you'd rather take the zero result and press on, but
I'd argue that they're subtle ones that you're not likely gonna be writing
in SQL.

Anyway, for now I propose the attached patch.  The test cases inquire into
the edge-case behavior of pow() much more closely than we have done in the
past, and I wouldn't be a bit surprised if some of the older buildfarm
critters fail some of them.  So my inclination is to try this only in
HEAD for starters.  Even if we want to back-patch, I'd be hesitant to
put it in versions older than v12, where we started to require C99.

            regards, tom lane

diff --git a/src/backend/utils/adt/float.c b/src/backend/utils/adt/float.c
index 6a717f19bb..84d37de930 100644
--- a/src/backend/utils/adt/float.c
+++ b/src/backend/utils/adt/float.c
@@ -1565,7 +1565,7 @@ dpow(PG_FUNCTION_ARGS)

     if (unlikely(isinf(result)) && !isinf(arg1) && !isinf(arg2))
         float_overflow_error();
-    if (unlikely(result == 0.0) && arg1 != 0.0)
+    if (unlikely(result == 0.0) && arg1 != 0.0 && !isinf(arg1) && !isinf(arg2))
         float_underflow_error();

     PG_RETURN_FLOAT8(result);
@@ -1581,15 +1581,38 @@ dexp(PG_FUNCTION_ARGS)
     float8        arg1 = PG_GETARG_FLOAT8(0);
     float8        result;

-    errno = 0;
-    result = exp(arg1);
-    if (errno == ERANGE && result != 0 && !isinf(result))
-        result = get_float8_infinity();
-
-    if (unlikely(isinf(result)) && !isinf(arg1))
-        float_overflow_error();
-    if (unlikely(result == 0.0))
-        float_underflow_error();
+    /*
+     * Handle NaN and Inf cases explicitly.  This avoids needing to assume
+     * that the platform's exp() conforms to POSIX for these cases, and it
+     * removes some edge cases for the overflow checks below.
+     */
+    if (isnan(arg1))
+        result = arg1;
+    else if (isinf(arg1))
+    {
+        /* Per POSIX, exp(-Inf) is 0 */
+        result = (arg1 > 0.0) ? arg1 : 0;
+    }
+    else
+    {
+        /*
+         * On some platforms, exp() will not set errno but just return Inf or
+         * zero to report overflow/underflow; therefore, test both cases.
+         */
+        errno = 0;
+        result = exp(arg1);
+        if (unlikely(errno == ERANGE))
+        {
+            if (result != 0.0)
+                float_overflow_error();
+            else
+                float_underflow_error();
+        }
+        else if (unlikely(isinf(result)))
+            float_overflow_error();
+        else if (unlikely(result == 0.0))
+            float_underflow_error();
+    }

     PG_RETURN_FLOAT8(result);
 }
diff --git a/src/test/regress/expected/float8.out b/src/test/regress/expected/float8.out
index aaef20bcfd..3957fb58d8 100644
--- a/src/test/regress/expected/float8.out
+++ b/src/test/regress/expected/float8.out
@@ -385,6 +385,158 @@ SELECT power(float8 'NaN', float8 '0');
      1
 (1 row)

+SELECT power(float8 'inf', float8 '0');
+ power
+-------
+     1
+(1 row)
+
+SELECT power(float8 '-inf', float8 '0');
+ power
+-------
+     1
+(1 row)
+
+SELECT power(float8 '0', float8 'inf');
+ power
+-------
+     0
+(1 row)
+
+SELECT power(float8 '0', float8 '-inf');
+ERROR:  zero raised to a negative power is undefined
+SELECT power(float8 '1', float8 'inf');
+ power
+-------
+     1
+(1 row)
+
+SELECT power(float8 '1', float8 '-inf');
+ power
+-------
+     1
+(1 row)
+
+SELECT power(float8 '-1', float8 'inf');
+ power
+-------
+     1
+(1 row)
+
+SELECT power(float8 '-1', float8 '-inf');
+ power
+-------
+     1
+(1 row)
+
+SELECT power(float8 '0.1', float8 'inf');
+ power
+-------
+     0
+(1 row)
+
+SELECT power(float8 '-0.1', float8 'inf');
+ power
+-------
+     0
+(1 row)
+
+SELECT power(float8 '1.1', float8 'inf');
+  power
+----------
+ Infinity
+(1 row)
+
+SELECT power(float8 '-1.1', float8 'inf');
+  power
+----------
+ Infinity
+(1 row)
+
+SELECT power(float8 '0.1', float8 '-inf');
+  power
+----------
+ Infinity
+(1 row)
+
+SELECT power(float8 '-0.1', float8 '-inf');
+  power
+----------
+ Infinity
+(1 row)
+
+SELECT power(float8 '1.1', float8 '-inf');
+ power
+-------
+     0
+(1 row)
+
+SELECT power(float8 '-1.1', float8 '-inf');
+ power
+-------
+     0
+(1 row)
+
+SELECT power(float8 'inf', float8 '-2');
+ power
+-------
+     0
+(1 row)
+
+SELECT power(float8 'inf', float8 '2');
+  power
+----------
+ Infinity
+(1 row)
+
+SELECT power(float8 'inf', float8 'inf');
+  power
+----------
+ Infinity
+(1 row)
+
+SELECT power(float8 'inf', float8 '-inf');
+ power
+-------
+     0
+(1 row)
+
+SELECT power(float8 '-inf', float8 '-2');
+ power
+-------
+     0
+(1 row)
+
+SELECT power(float8 '-inf', float8 '-3');
+ power
+-------
+    -0
+(1 row)
+
+SELECT power(float8 '-inf', float8 '2');
+  power
+----------
+ Infinity
+(1 row)
+
+SELECT power(float8 '-inf', float8 '3');
+   power
+-----------
+ -Infinity
+(1 row)
+
+SELECT power(float8 '-inf', float8 'inf');
+  power
+----------
+ Infinity
+(1 row)
+
+SELECT power(float8 '-inf', float8 '-inf');
+ power
+-------
+     0
+(1 row)
+
 -- take exp of ln(f.f1)
 SELECT '' AS three, f.f1, exp(ln(f.f1)) AS exp_ln_f1
    FROM FLOAT8_TBL f
@@ -396,6 +548,13 @@ SELECT '' AS three, f.f1, exp(ln(f.f1)) AS exp_ln_f1
        | 1.2345678901234e-200 | 1.23456789012339e-200
 (3 rows)

+-- check edge cases for exp
+SELECT exp('inf'::float8), exp('-inf'::float8), exp('nan'::float8);
+   exp    | exp | exp
+----------+-----+-----
+ Infinity |   0 | NaN
+(1 row)
+
 -- cube root
 SELECT ||/ float8 '27' AS three;
  three
diff --git a/src/test/regress/sql/float8.sql b/src/test/regress/sql/float8.sql
index e540f03b07..3a8c737fb2 100644
--- a/src/test/regress/sql/float8.sql
+++ b/src/test/regress/sql/float8.sql
@@ -120,12 +120,41 @@ SELECT power(float8 'NaN', float8 'NaN');
 SELECT power(float8 '-1', float8 'NaN');
 SELECT power(float8 '1', float8 'NaN');
 SELECT power(float8 'NaN', float8 '0');
+SELECT power(float8 'inf', float8 '0');
+SELECT power(float8 '-inf', float8 '0');
+SELECT power(float8 '0', float8 'inf');
+SELECT power(float8 '0', float8 '-inf');
+SELECT power(float8 '1', float8 'inf');
+SELECT power(float8 '1', float8 '-inf');
+SELECT power(float8 '-1', float8 'inf');
+SELECT power(float8 '-1', float8 '-inf');
+SELECT power(float8 '0.1', float8 'inf');
+SELECT power(float8 '-0.1', float8 'inf');
+SELECT power(float8 '1.1', float8 'inf');
+SELECT power(float8 '-1.1', float8 'inf');
+SELECT power(float8 '0.1', float8 '-inf');
+SELECT power(float8 '-0.1', float8 '-inf');
+SELECT power(float8 '1.1', float8 '-inf');
+SELECT power(float8 '-1.1', float8 '-inf');
+SELECT power(float8 'inf', float8 '-2');
+SELECT power(float8 'inf', float8 '2');
+SELECT power(float8 'inf', float8 'inf');
+SELECT power(float8 'inf', float8 '-inf');
+SELECT power(float8 '-inf', float8 '-2');
+SELECT power(float8 '-inf', float8 '-3');
+SELECT power(float8 '-inf', float8 '2');
+SELECT power(float8 '-inf', float8 '3');
+SELECT power(float8 '-inf', float8 'inf');
+SELECT power(float8 '-inf', float8 '-inf');

 -- take exp of ln(f.f1)
 SELECT '' AS three, f.f1, exp(ln(f.f1)) AS exp_ln_f1
    FROM FLOAT8_TBL f
    WHERE f.f1 > '0.0';

+-- check edge cases for exp
+SELECT exp('inf'::float8), exp('-inf'::float8), exp('nan'::float8);
+
 -- cube root
 SELECT ||/ float8 '27' AS three;