Closed Bug 1775254 Opened 2 years ago Closed 2 years ago

Improve Math.pow accuracy for large exponents

Categories

(Core :: JavaScript Engine, enhancement, P3)

Firefox 101
enhancement

Tracking

()

RESOLVED FIXED
104 Branch
Tracking Status
firefox104 --- fixed

People

(Reporter: maxgraey, Assigned: anba)

References

(Blocks 1 open bug)

Details

Attachments

(2 files)

User Agent: Mozilla/5.0 (Macintosh; Intel Mac OS X 10.15; rv:101.0) Gecko/20100101 Firefox/101.0

Steps to reproduce:

For power with exponent >= 64 we got pretty imprecise results

Actual results:

Math.pow(10, 63)

1e+63 (ok)

Math.pow(10, 64)

1.0000000000000002e+64

Math.pow(10, 208)

1.0000000000000004e+208

Math.pow(10, 308)

1.0000000000000006e+308

Expected results:

Math.pow(10, 64)

1e+64

Math.pow(10, 208)

1e+208

Math.pow(10, 308)

1e+308

These are the results we can get in Chrome

Moving to a more appropriate component.

Component: Untriaged → JavaScript Engine
Product: Firefox → Core

I strongly suspect this isn't an accuracy issue, but rather a printing issue (which is to say, Chrome is printing the same value differently).

Safari agrees with Firefox on how it prints these.

Note as well, the JavaScript Number type is only accurate up to Number.MAX_SAFE_INTEGER, and after that you're going to have to be aware of floating point inaccuracy problems.

If working with integers of this size, consider BigInt values:

10n ** 64n

100000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000n

Severity: -- → N/A
Status: UNCONFIRMED → RESOLVED
Closed: 2 years ago
Priority: -- → P5
Resolution: --- → INVALID

I strongly suspect this isn't an accuracy issue, but rather a printing issue (which is to say, Chrome is printing the same value differently).

That's not true. Just to check this in Firefox:

Math.pow(10, 208) == 1e208
false

and in Chrome:

Math.pow(10, 208) == 1e208
true

Severity: N/A → --
Status: RESOLVED → REOPENED
Ever confirmed: true
Priority: P5 → P3
Resolution: INVALID → ---

So, on this test case:

for (var i = 0; i < 200; i++) { 
  var res = eval(`Math.pow(10,${i}) == 1e${i}`);
  if (!res) {
    console.log(`${i}: ${Math.pow(10,i)} != ${eval("1e"+i)}`)
  }
} 

Everything goes sideways in both firefox and safari at i==33.

However: similar issues exist for Chrome (Canary) too:

VM10:3 26: 9.999999999999999e+25 != 1e+26
VM10:3 29: 1.0000000000000001e+29 != 1e+29
VM10:3 34: 1.0000000000000001e+34 != 1e+34
VM10:3 36: 9.999999999999999e+35 != 1e+36
VM10:3 39: 1.0000000000000001e+39 != 1e+39
VM10:3 42: 9.999999999999999e+41 != 1e+42
VM10:3 49: 1.0000000000000001e+49 != 1e+49
VM10:3 58: 1.0000000000000001e+58 != 1e+58
VM10:3 61: 1.0000000000000001e+61 != 1e+61
VM10:3 70: 9.999999999999999e+69 != 1e+70
VM10:3 95: 9.999999999999999e+94 != 1e+95
VM10:3 96: 9.999999999999999e+95 != 1e+96
VM10:3 111: 1.0000000000000001e+111 != 1e+111
VM10:3 126: 1.0000000000000001e+126 != 1e+126
VM10:3 136: 9.999999999999999e+135 != 1e+136
VM10:3 142: 9.999999999999999e+141 != 1e+142
VM10:3 165: 1.0000000000000001e+165 != 1e+165
VM10:3 166: 1.0000000000000001e+166 != 1e+166
VM10:3 169: 1.0000000000000001e+169 != 1e+169
VM10:3 187: 1.0000000000000001e+187 != 1e+187
VM10:3 191: 9.999999999999999e+190 != 1e+191

This is outside my numerical methods expertise, but sure looks like alternative algorithms producing differing roundings.

I'm not sure the spec ever guarantees that you'd get 1eX == Math.pow(10,X)

Status: REOPENED → RESOLVED
Closed: 2 years ago2 years ago
Resolution: --- → WONTFIX

Informal all maths libs should have guarantee ULP == 1 or at least ULP < 1.5 for all functions. But Math.pow with exp > 64 breaks this rule significantly. For example:

function toU64(x) {
const
U64 = new BigUint64Array(1),
F64 = new Float64Array(U64.buffer);
F64[0] = x;
return U64[0];
}

function absU64(x) {
return x < 0n
? BigInt.asUintN(64, -x)
: x;
}

function ulpDistance(x, y) {
const
ux = toU64(x),
uy = toU64(y);
return Number(absU64(absU64(ux) - absU64(uy)));
}

const actual = Math.pow(10, 208); // 1.0000000000000004e+208 on Firefox
const ideal = 1e208;
console.log(ulpDistance(actual, ideal) + ' ulp');

it returns 4 ulp on Firefox!

Please note Math.pow(10, 208); and 1e208 were taken as an example to demonstrate a significant deviation from the desired accuracy. ULP > 2 means that the mathematical function is not suitable for scientific purposes because its accuracy is worse than the rounding error, which is considered a bad sign. I suggest that you check the accuracy of the same function in industrial languages such as C++, Rust, Python or Julia. Chrome has a more precise implementation of the power function. Why not just follow their example?

By the spec, what we have is fine. Math.pow is specified to "Return an implementation-approximated Number value representing the result of raising ℝ(base) to the ℝ(exponent) power."

I agree, it would be nice to improve numerical accuracy here; however, given this isn't incorrect by spec, marking as a defect an enhancement.

Severity: -- → N/A
Status: RESOLVED → REOPENED
Type: defect → enhancement
Resolution: WONTFIX → ---
Summary: Inaccurate Math.pow for huge arguments → Improve Math.pow accuracy for large exponents

(I could be wrong about SpiderMonkey using fdlibm::pow -- at least once upon a time the idea was that fdlibm calls were namespace prefixed; however, I did try changing that pow call to use fdlibm::pow explicitly, and got the same ulp errors)

Hmm, but Rust also uses libm implementation (https://github.com/rust-lang/compiler-builtins). But it outputs correct result:

println!("{:e}", 10_f64.powf(208_f64));
1e208

I checked in others languages and all the same. It looks like FF uses own or outdated version of libm

Btw I don't recommend to use fdlibm. It was designed for small size but not for fast and precise implementation. Also I found similar issue which addressed to fdlibm in general: https://bugzilla.mozilla.org/show_bug.cgi?id=933257. But it closed

For anyone doing archaeology: Previous discussions

Interestingly, the specification actually recommends using fdlibm (oh, but which of the forks scattered to the wind?)

(In reply to Matthew Gaudet (he/him) [:mgaudet] from comment #8)

(I could be wrong about SpiderMonkey using fdlibm::pow -- at least once upon a time the idea was that fdlibm calls were namespace prefixed; however, I did try changing that pow call to use fdlibm::pow explicitly, and got the same ulp errors)

Yes, we don't use fdlibm::pow because it is slower than std::pow (see bug 933257). std::pow from libc has maximal error of 1 ULP [1], whereas fdlibm::pow returns the "nearly rounded" result [2]. 1 ULP error was considered acceptable, therefore we went with using std::pow.

Anyway, the issue from this bug isn't actually related to std::pow, but instead powi. The error from powi can be quite large, for example try (-1.00000041289256280663266807096078991889953613281250000) ** -365287834. JSC shows the same behaviour, because they have a similar optimisation [3].

[1] https://www.gnu.org/software/libc/manual/html_node/Errors-in-Math-Functions.html#Errors-in-Math-Functions
[2] https://github.com/freebsd/freebsd-src/blob/master/lib/msun/src/e_pow.c
[3] https://github.com/WebKit/WebKit/blob/3a8b3ba34d64ecbb737ee6c4b0a7fd86c0147452/Source/JavaScriptCore/runtime/MathCommon.cpp#L436-L448

I see, thanks!
v8 recently introduced another wrong optimization:
https://bugs.chromium.org/p/v8/issues/detail?id=12996&q=pow&can=2

In the pursuit of good numbers in benchmarks every browser systematically deteriorates the accuracy of mathematics, which can subsequently lead to disastrous consequences. I think that TC39 we should declare error bounds for math functions at 1.5 ULP level in spec and add according tests to check this. That nobody would be tempted to make optimizations worsening accuracy.

The previous js::powi implementation could cause large floating point
precision errors. For example (-1.00000041289256280663266807096078991889953613281250000) ** -365287834
returned 3.1456398721089196e-66 instead of 3.1456398910784315e-66. This is an
error of 35987774 ULP.

Change js::powi as follows:

  • Optimise only for non-negative exponents.
  • Always optimise when the exponent is <= 4 to match MPow::foldsTo.
  • When the exponent is > 4, only optimise when the base is an integer. Perform
    the whole computation with integer values and when an overflow was detected,
    fallback to std::pow.
  • In all other cases use std::pow.
Assignee: nobody → andrebargull

"jsmath.cpp" includes <cmath>, so we should include the "std" namespace.

Depends on D150791

Change js::powi as follows:

Optimise only for non-negative exponents.
Always optimise when the exponent is <= 4 to match MPow::foldsTo.
When the exponent is > 4, only optimise when the base is an integer. Perform
the whole computation with integer values and when an overflow was detected,
fallback to std::pow.
In all other cases use std::pow.

x * x * x still have at least 2 Ulp error compare to x ** 3, x * x * x * x even more.

For Chrome's team I proposed trade-off approach described in "Computing Integer Powers in Floating-Point Arithmetic" paper:
https://arxiv.org/pdf/0705.4369.pdf

This should maintain accuracy and at the same time be much faster than a powf(x, 3.0)

Here are working proof of concept for fast x^3 on Rust:
https://play.rust-lang.org/?version=nightly&mode=debug&edition=2021&gist=674dc6e818fa0fa7122328d24e62cf3c

WDYT?

(In reply to MaxGraey from comment #17)

For Chrome's team I proposed trade-off approach described in "Computing Integer Powers in Floating-Point Arithmetic" paper:
https://arxiv.org/pdf/0705.4369.pdf

The algorithms in the paper assume that no subnormal numbers are used. From the "Introduction" paragraph:

[...] and that the input as well as the output values of the power function are not subnormal numbers, [...]

So I don't think we can use the algorithms proposed in that paper, right?

Apart from that, std::pow can also have an error larger than 2 ULP on some platforms (= Windows) [1], so I'm not too worried about the x ** 3 case.


[1] -0.99999908129426284819629699995857663452625274658203125 ** -172780371 returns -8.659047226014985e+68 (0xce400f20e3e56451) on Windows, whereas fdlibm and musl both compute -8.65904722601499e+68 (0xce400f20e3e56454).

So instead of fixing the situation in Windows and making JavaScript deterministic at least on a cross-platform level and using the power function from musl for example, you are just proposing to degrade accuracy elsewhere as well? This seems to be the Broken Windows Theory (https://en.wikipedia.org/wiki/Broken_windows_theory) in action, don't you think?

What I see as the main problem with speculative optimizations like x ** 3 -> x * x * x by compiler. Users who care about performance are more likely to explicitly do it themselves. But users who care about accuracy now have no easy way to improve their code, because now they can't rely on Math.pow in this case.

The algorithms in the paper assume that no subnormal numbers are used. From the "Introduction" paragraph:

I don't think it's a problem to insert a check like:

if (isSubnormal(x)) {
// use fast path
} else {
// fallback to full precision power
}

Or you could put the fastpath directly into the pow function, but it would have to be completely in-house (e.g., based on the musl implementation)

It is important from a fingerprinting angle to respect current usage of fdlibm to ensure consistent behavior across different OSes/versions/hardware. The recently-ish added pref javascript.options.use_fdlibm_for_sin_cos_tan improved our behavior and I don't think this patch undoes that, but I also want to ensure we don't go backwards on cross-platform consistency for other functions like pow.

(In reply to Tom Ritter [:tjr] from comment #22)

[...] but I also want to ensure we don't go backwards on cross-platform consistency for other functions like pow.

We were already calling std::pow in other cases, with this patch it's only called more often. Probably std::pow was overlooked in bug 531915?

Probably std::pow was overlooked in bug 531915?

OT: That patch fixed AFAICT desktop entropy, and reduced it in Android, but yes there are others. Bug 1722181 mentions links to pow, but it's not the only one - this is Android: I do not know of any other math entropy on desktop with the RFP fdlibm sin/cos/tan patch

Pushed by andre.bargull@gmail.com:
https://hg.mozilla.org/integration/autoland/rev/635668eb1697
Part 1: Avoid large floating point precision errors in js::powi. r=mgaudet
https://hg.mozilla.org/integration/autoland/rev/546f112a4591
Part 2: Prefix cmath calls with "std::". r=mgaudet
Status: REOPENED → RESOLVED
Closed: 2 years ago2 years ago
Resolution: --- → FIXED
Target Milestone: --- → 104 Branch
Regressions: 1784111
See Also: → 1798789
You need to log in before you can comment on or make changes to this bug.

Attachment

General

Created:
Updated:
Size: