Exponentiation confusion

On 13/10/2022 18:20 CEST Adrian Klaver <adrian.klaver@aklaver.com> wrote:

In trying to answer an SO question I ran across this:

Postgres version 14.5

select 10^(-1 * 18);
?column?
----------
1e-18

select 10^(-1 * 18::numeric);
?column?
--------------------
0.0000000000000000

Same for power:

select power(10, -18);
power
-------
1e-18
(1 row)

select power(10, -18::numeric);
power
--------------------
0.0000000000000000

Why is the cast throwing off the result?

power has two overloads: https://www.postgresql.org/docs/14/functions-math.html#id-1.5.8.9.6.2.2.19.1.1.1

Calling power(numeric, numeric) is what I expect in that case instead of
downcasting the exponent argument to double precision, thus losing precision.

select
pg_typeof(power(10, -18)),
pg_typeof(power(10, -18::numeric));

pg_typeof | pg_typeof
------------------+-----------
double precision | numeric
(1 row)

Determining the right function is described in https://www.postgresql.org/docs/14/typeconv-func.html

--
Erik

On 2022-10-13 09:20:51 -0700, Adrian Klaver wrote:

In trying to answer an SO question I ran across this:

Postgres version 14.5

Same for 11.17. So it's been like that for some time, maybe forever.

select power(10, -18);
power
-------
1e-18
(1 row)

select power(10, -18::numeric);
power
--------------------
0.0000000000000000

Why is the cast throwing off the result?

It seems that the number of decimals depends only on the first argument:

hjp=> select power(10::numeric, -2::numeric);
╔════════════════════╗
║ power ║
╟────────────────────╢
║ 0.0100000000000000 ║
╚════════════════════╝
(1 row)
hjp=> select power(10::numeric, -16::numeric);
╔════════════════════╗
║ power ║
╟────────────────────╢
║ 0.0000000000000001 ║
╚════════════════════╝
(1 row)
hjp=> select power(10::numeric, -18::numeric);
╔════════════════════╗
║ power ║
╟────────────────────╢
║ 0.0000000000000000 ║
╚════════════════════╝
(1 row)

hjp=> select power(10::numeric, 18::numeric);
╔══════════════════════════════════════╗
║ power ║
╟──────────────────────────────────────╢
║ 1000000000000000000.0000000000000000 ║
╚══════════════════════════════════════╝
(1 row)

hjp=> select power(10::numeric(32,30), 18::numeric);
╔════════════════════════════════════════════════════╗
║ power ║
╟────────────────────────────────────────────────────╢
║ 1000000000000000000.000000000000000000000000000000 ║
╚════════════════════════════════════════════════════╝
(1 row)
hjp=> select power(10::numeric(32,30), -16::numeric);
╔══════════════════════════════════╗
║ power ║
╟──────────────────────────────────╢
║ 0.000000000000000100000000000000 ║
╚══════════════════════════════════╝
(1 row)

So the number of decimals by default isn't sufficient to represent
10^-18. You have to explicitely increase it.

hp

--
_ | Peter J. Holzer | Story must make more sense than reality.
|_|_) | |
| | | hjp@hjp.at | -- Charles Stross, "Creative writing
__/ | http://www.hjp.at/ | challenge!"

Erik Wienhold <ewie@ewie.name> writes:

On 13/10/2022 18:20 CEST Adrian Klaver <adrian.klaver@aklaver.com> wrote:

select power(10, -18::numeric);
power
--------------------
0.0000000000000000

Why is the cast throwing off the result?

Calling power(numeric, numeric) is what I expect in that case instead of
downcasting the exponent argument to double precision, thus losing precision.

An inexact result isn't surprising, but it shouldn't be *that* inexact.
It looks to me like numeric.c's power_var_int() code path is setting the
result rscale without considering the possibility that the result will
have negative weight (i.e. be less than one). The main code path in
power_var() does adjust for that, so for example

regression=# select power(10, -18.00000001::numeric);
power
-------------------------------------
0.000000000000000000999999976974149
(1 row)

but with an exact-integer exponent, not so much --- you just get 16 digits
which isn't enough.

I'm inclined to think that we should push the responsibility for choosing
its rscale into power_var_int(), because internally that already does
estimate the result weight, so with a little code re-ordering we won't
need duplicative estimates. Don't have time to work on that right now
though ... Dean, are you interested in fixing this?

regards, tom lane

On Thu, 13 Oct 2022 at 18:17, Tom Lane <tgl@sss.pgh.pa.us> wrote:

I'm inclined to think that we should push the responsibility for choosing
its rscale into power_var_int(), because internally that already does
estimate the result weight, so with a little code re-ordering we won't
need duplicative estimates. Don't have time to work on that right now
though ... Dean, are you interested in fixing this?

OK, I'll take a look.

The most obvious thing to do is to try to make power_var_int() choose
the same result rscale as power_var() so that the results are
consistent regardless of whether the exponent is an integer.

It's worth noting, however, that that will cause in a *reduction* in
the output rscale rather than an increase in some cases, since the
power_var_int() code path currently always chooses an rscale of at
least 16, whereas the other code path in power_var() uses the rscales
of the 2 inputs, and produces a minimum of 16 significant digits,
rather than 16 digits after the decimal point. For example:

select power(5.678, 18.00000001::numeric);
power
-------------------------
37628507689498.14987457
(1 row)

select power(5.678, 18::numeric);
power
---------------------------------
37628507036041.8454541428979479
(1 row)

Regards,
Dean

Dean Rasheed <dean.a.rasheed@gmail.com> writes:

The most obvious thing to do is to try to make power_var_int() choose
the same result rscale as power_var() so that the results are
consistent regardless of whether the exponent is an integer.

Yeah, I think we should try to end up with that.

It's worth noting, however, that that will cause in a *reduction* in
the output rscale rather than an increase in some cases, since the
power_var_int() code path currently always chooses an rscale of at
least 16, whereas the other code path in power_var() uses the rscales
of the 2 inputs, and produces a minimum of 16 significant digits,
rather than 16 digits after the decimal point.

Right. I think this is not bad though. In a lot of cases (such
as the example here) the current behavior is just plastering on
useless zeroes.

regards, tom lane

On 13/10/2022 19:16 CEST Tom Lane <tgl@sss.pgh.pa.us> wrote:

Erik Wienhold <ewie@ewie.name> writes:

On 13/10/2022 18:20 CEST Adrian Klaver <adrian.klaver@aklaver.com> wrote:

select power(10, -18::numeric);
power
--------------------
0.0000000000000000

Why is the cast throwing off the result?

Calling power(numeric, numeric) is what I expect in that case instead of
downcasting the exponent argument to double precision, thus losing precision.

An inexact result isn't surprising, but it shouldn't be *that* inexact.

Ah, now I see the problem. I saw a bunch of zeros but not that it's *all*
zeros. Nevermind.

--
Erik

[Moving this to -hackers]

On 13/10/2022 18:20 CEST Adrian Klaver <adrian.klaver@aklaver.com> wrote:

select power(10, -18::numeric);
power
--------------------
0.0000000000000000

On Thu, 13 Oct 2022 at 18:17, Tom Lane <tgl@sss.pgh.pa.us> wrote:

An inexact result isn't surprising, but it shouldn't be *that* inexact.

I'm inclined to think that we should push the responsibility for choosing
its rscale into power_var_int(), because internally that already does
estimate the result weight, so with a little code re-ordering we won't
need duplicative estimates. Don't have time to work on that right now
though ... Dean, are you interested in fixing this?

Here's a patch along those lines, bringing power_var_int() more in
line with neighbouring functions by having it choose its own result
scale.

It was necessary to also move the overflow/underflow tests up, in
order to avoid a potential integer overflow when deciding the rscale.

Looking more closely at the upper limit of the overflow test, it turns
out it was far too large. I'm not sure where the "3 * SHRT_MAX" came
from, but I suspect it was just a thinko on my part, back in
7d9a4737c2. I've replaced that with SHRT_MAX + 1, which kicks in much
sooner without changing the actual maximum result allowed, which is <
10^131072 (the absolute upper limit of the numeric type).

The first half the the underflow test condition "f + 1 < -rscale" goes
away, since this is now being done before rscale is computed, and the
choice of rscale makes that condition false. In fact, the new choice
of rscale now ensures that when sig_digits is computed further down,
it is guaranteed to be strictly greater than 0, rather than merely
being >= 0 as before, which is good.

As expected, various regression test results change, since the number
of significant digits computed is now different, but I think the new
results look a lot better, and more consistent. I regenerated the
numeric_big test results by re-running the bc script and rounding to
the new output precisions, and the results from power_var_int()
exactly match in every case. This already included a number of cases
that used to round to zero, and now produce much more reasonable
results.

The test cases where the result actually does round to zero now output
1000 zeros after the decimal point. That looks a little messy, but I
think it's the right thing to do in fixed-point arithmetic -- it's
consistent with the fractional power case, and with exp(numeric),
reflecting the fact that the result is zero to 1000 decimal places,
whilst not being exactly zero.

Overall, I'm quite happy with these results. The question is, should
this be back-patched?

In the past, I think I've only back-patched numeric bug-fixes where
the digits output by the old code were incorrect or an error was
thrown, not changes that resulted in a different number of digits
being output, changing the precision of already-correct results.
However, having 10.0^(-18) produce zero seems pretty bad, so my
inclination is to back-patch, unless anyone objects.

Regards,
Dean

On Tue, Oct 18, 2022 at 6:18 AM Dean Rasheed <dean.a.rasheed@gmail.com> wrote:

Overall, I'm quite happy with these results. The question is, should
this be back-patched?

In the past, I think I've only back-patched numeric bug-fixes where
the digits output by the old code were incorrect or an error was
thrown, not changes that resulted in a different number of digits
being output, changing the precision of already-correct results.
However, having 10.0^(-18) produce zero seems pretty bad, so my
inclination is to back-patch, unless anyone objects.

I don't think that back-patching is a very good idea. The bar for
changing query results should be super-high. Applications can depend
on the existing behavior even if it's wrong.

--
Robert Haas
EDB: http://www.enterprisedb.com

On Tue, 18 Oct 2022 at 20:18, Robert Haas <robertmhaas@gmail.com> wrote:

On Tue, Oct 18, 2022 at 6:18 AM Dean Rasheed <dean.a.rasheed@gmail.com> wrote:

Overall, I'm quite happy with these results. The question is, should
this be back-patched?

In the past, I think I've only back-patched numeric bug-fixes where
the digits output by the old code were incorrect or an error was
thrown, not changes that resulted in a different number of digits
being output, changing the precision of already-correct results.
However, having 10.0^(-18) produce zero seems pretty bad, so my
inclination is to back-patch, unless anyone objects.

I don't think that back-patching is a very good idea. The bar for
changing query results should be super-high. Applications can depend
on the existing behavior even if it's wrong.

OK, on reflection, I think that makes sense. Applied to HEAD only.

Regards,
Dean