julia icon indicating copy to clipboard operation
julia copied to clipboard

speed up Ryu.pow5

Open oscardssmith opened this issue 2 years ago • 3 comments

The current algorithm does up to p divisions, while this algorithm does only 1. I've assigned @quinnj to review since I'm not 100% sure that there isn't a very good reason that this function isn't written like this.

oscardssmith avatar Sep 14 '22 16:09 oscardssmith

I'm not 100% sure that there isn't a very good reason that this function isn't written like this.

My guess would be overflow protection for sufficiently large p? Does wraparound prevent any shenanigans that could happen due to that?

Seelengrab avatar Sep 15 '22 06:09 Seelengrab

Do you know if this needs an overflow check?

oscardssmith avatar Sep 15 '22 18:09 oscardssmith

I don't off the top of my head

quinnj avatar Sep 15 '22 18:09 quinnj

Maybe run a PkgEval and merge if it looks ok?

KristofferC avatar Sep 22 '22 11:09 KristofferC

It feels kind of silly to use pkgeval for something this small, but I guess we might as well.

oscardssmith avatar Sep 22 '22 14:09 oscardssmith

Or just merge, tests pass after all and I guess we will find it in release PkgEval if it is too bad

KristofferC avatar Sep 22 '22 14:09 KristofferC

#yolo

quinnj avatar Oct 06 '22 06:10 quinnj

Felt I should look at this since printing numbers incorrectly would be bad. This use of pow5 seems safe:

https://github.com/JuliaLang/julia/blob/f71b8398e6f6eca93eb6f7cefaa940aa7e2b1f5a/base/ryu/shortest.jl#L58-L66

My reasoning is that 5^p is correct for p ≤ 27 and here we have that q ≤ qinvbound(T) defined here:

https://github.com/JuliaLang/julia/blob/f71b8398e6f6eca93eb6f7cefaa940aa7e2b1f5a/base/ryu/utils.jl#L17-L19

You can see that q can be at most 21, so it's safe. The other use of pow5 is less clear:

https://github.com/JuliaLang/julia/blob/f71b8398e6f6eca93eb6f7cefaa940aa7e2b1f5a/base/ryu/exp.jl#L152-L159

Does anyone know what the range of possible values for precision and e are here? If e can be more than 27 greater than precision then we might have a problem with this change.

StefanKarpinski avatar Oct 27 '22 17:10 StefanKarpinski

Bump: @oscardssmith, @quinnj do you know about the range of values that e can take here?

StefanKarpinski avatar Nov 08 '22 15:11 StefanKarpinski

ah it appears that e can possibly be large if someone writes a decimal with a ton of digits.

oscardssmith avatar Nov 08 '22 16:11 oscardssmith

Do you have an example?

KristofferC avatar Nov 08 '22 17:11 KristofferC

julia> using Printf

julia> @printf "%.8g" 4.645833859177319e63 # master
4.6458338e+63

julia> @printf "%.8g" 4.645833859177319e63 # 1.8
4.6458339e+63

There are two ways that pow5 can be wrong,

  1. returning false when m2 is divisible by big(5)^requiredFives
  2. returning true when m2 is not divisible by big(5)^requiredFives In the first case, for pow5 to be wrong, big(5)^requiredFives must be greater than typemax(Int), but m2 is always less than 1<<53, and a small number is not divisible by a large number. Thus m2 will never be divisible by the true big(5)^requiredFives when there is overflow and case 1 is okay.

For the second case, we need to find a requiredFives such that pow5 may wrongly return true. That is, find a requiredFives such that there exists m2 where m2 % 5^requiredFives == 0 For this to be the case, either m2 must be 0 (handled much earlier as a special case) or abs(5^requiredFives) <= m2 < 1<<53. Searching exhaustively from the cases that have overflow, findfirst(requiredFives -> abs(5^requiredFives) < 1<<53, 28:10^5)+27 == 55. Note that 5^55 < 0. If, on the other hand, we used unsinged(5)^requiredFives, then we'd have findfirst(requiredFives -> unsigned(5)^requiredFives < 1<<53, 28:10^5)+27 == 1048 and IIUC requiredFives is at most log10(floatmax(Float64)) == 308.25471555991675.

So using unsigned(5) should fix this and I benchmark it as comparable to signed exponentiation.

LilithHafner avatar Nov 09 '22 13:11 LilithHafner