Skip to content

Commit

Permalink
Document internal Ryu functions, remove duplicate functionality. (Jul…
Browse files Browse the repository at this point in the history
…iaLang#33406)

Also moved some functions from Math to Base.
  • Loading branch information
simonbyrne committed Sep 29, 2019
1 parent 0b55ce1 commit d036575
Show file tree
Hide file tree
Showing 5 changed files with 96 additions and 44 deletions.
5 changes: 5 additions & 0 deletions base/float.jl
Original file line number Diff line number Diff line change
Expand Up @@ -876,6 +876,11 @@ significand_mask(::Type{Float16}) = 0x03ff
for T in (Float16, Float32, Float64)
@eval significand_bits(::Type{$T}) = $(trailing_ones(significand_mask(T)))
@eval exponent_bits(::Type{$T}) = $(sizeof(T)*8 - significand_bits(T) - 1)
@eval exponent_bias(::Type{$T}) = $(Int(exponent_one(T) >> significand_bits(T)))
# maximum float exponent
@eval exponent_max(::Type{$T}) = $(Int(exponent_mask(T) >> significand_bits(T)) - exponent_bias(T))
# maximum float exponent without bias
@eval exponent_raw_max(::Type{$T}) = $(Int(exponent_mask(T) >> significand_bits(T)))
end

# integer size of float
Expand Down
11 changes: 2 additions & 9 deletions base/math.jl
Original file line number Diff line number Diff line change
Expand Up @@ -22,7 +22,8 @@ import .Base: log, exp, sin, cos, tan, sinh, cosh, tanh, asin,

using .Base: sign_mask, exponent_mask, exponent_one,
exponent_half, uinttype, significand_mask,
significand_bits, exponent_bits
significand_bits, exponent_bits, exponent_bias,
exponent_max, exponent_raw_max

using Core.Intrinsics: sqrt_llvm

Expand All @@ -38,14 +39,6 @@ end
"Complex(x)^y, or similar.")))
end

for T in (Float16, Float32, Float64)
@eval exponent_bias(::Type{$T}) = $(Int(exponent_one(T) >> significand_bits(T)))
# maximum float exponent
@eval exponent_max(::Type{$T}) = $(Int(exponent_mask(T) >> significand_bits(T)) - exponent_bias(T))
# maximum float exponent without bias
@eval exponent_raw_max(::Type{$T}) = $(Int(exponent_mask(T) >> significand_bits(T)))
end

# non-type specific math functions

"""
Expand Down
9 changes: 8 additions & 1 deletion base/ryu/Ryu.jl
Original file line number Diff line number Diff line change
@@ -1,10 +1,17 @@
module Ryu

import .Base: significand_bits, exponent_bits, exponent_bias

include("utils.jl")
include("shortest.jl")
include("fixed.jl")
include("exp.jl")

"""
Ryu.neededdigits(T)
Number of digits necessary to represent type `T` in fixed-precision decimal.
"""
neededdigits(::Type{Float64}) = 309 + 17
neededdigits(::Type{Float32}) = 39 + 9 + 2
neededdigits(::Type{Float16}) = 9 + 5 + 9
Expand Down Expand Up @@ -120,4 +127,4 @@ end

Base.print(io::IO, x::Union{Float16, Float32}) = show(io, x, true)

end # module
end # module
12 changes: 6 additions & 6 deletions base/ryu/shortest.jl
Original file line number Diff line number Diff line change
Expand Up @@ -75,15 +75,15 @@
return pos + neg + 3 + (typed && x isa Union{Float32, Float16} ? 2 : 0)
end

bits = uint(x)
mant = bits & (oftype(bits, 1) << mantissabits(T) - oftype(bits, 1))
exp = Int((bits >> mantissabits(T)) & ((Int64(1) << exponentbits(T)) - 1))
m2 = oftype(bits, Int64(1) << mantissabits(T)) | mant
e2 = exp - bias(T) - mantissabits(T)
bits = reinterpret(Unsigned, x)
mant = bits & (oftype(bits, 1) << significand_bits(T) - oftype(bits, 1))
exp = Int((bits >> significand_bits(T)) & ((Int64(1) << exponent_bits(T)) - 1))
m2 = oftype(bits, Int64(1) << significand_bits(T)) | mant
e2 = exp - exponent_bias(T) - significand_bits(T)
fraction = m2 & ((oftype(bits, 1) << -e2) - 1)
if e2 > 0 || e2 < -52 || fraction != 0
if exp == 0
e2 = 1 - bias(T) - mantissabits(T) - 2
e2 = 1 - exponent_bias(T) - significand_bits(T) - 2
m2 = mant
else
e2 -= 2
Expand Down
103 changes: 75 additions & 28 deletions base/ryu/utils.jl
Original file line number Diff line number Diff line change
@@ -1,25 +1,9 @@
const MANTISSA_MASK = 0x000fffffffffffff
const EXP_MASK = 0x00000000000007ff
const MANTISSA_MASK = Base.significand_mask(Float64)
const EXP_MASK = Base.exponent_mask(Float64) >> Base.significand_bits(Float64)

memcpy(d, doff, s, soff, n) = ccall(:memcpy, Cvoid, (Ptr{Cvoid}, Ptr{Cvoid}, Int), d + doff - 1, s + soff - 1, n)
memmove(d, doff, s, soff, n) = ccall(:memmove, Cvoid, (Ptr{Cvoid}, Ptr{Cvoid}, Int), d + doff - 1, s + soff - 1, n)

uint(x::Float16) = Core.bitcast(UInt16, x)
uint(x::Float32) = Core.bitcast(UInt32, x)
uint(x::Float64) = Core.bitcast(UInt64, x)

mantissabits(::Type{Float16}) = 10
mantissabits(::Type{Float32}) = 23
mantissabits(::Type{Float64}) = 52

exponentbits(::Type{Float16}) = 5
exponentbits(::Type{Float32}) = 8
exponentbits(::Type{Float64}) = 11

bias(::Type{Float16}) = 15
bias(::Type{Float32}) = 127
bias(::Type{Float64}) = 1023

pow5_bitcount(::Type{Float16}) = 30
pow5_bitcount(::Type{Float32}) = 61
pow5_bitcount(::Type{Float64}) = 121
Expand All @@ -36,16 +20,58 @@ qbound(::Type{Float16}) = 15
qbound(::Type{Float32}) = 31
qbound(::Type{Float64}) = 63

"""
Ryu.log10pow2(e::Integer)
Computes `floor(log10(2^e))`. This is valid for all `e < 1651`.
"""
log10pow2(e) = (e * 78913) >> 18


"""
Ryu.log10pow5(e::Integer)
Computes `floor(log10(5^e))`. This is valid for all `e < 2621`.
"""
log10pow5(e) = (e * 732923) >> 20

"""
Ryu.pow5bits(e)
Computes `e == 0 ? 1 : ceil(log2(5^e))`. This is valid for `e < 3529` (if performend in `Int32` arithmetic).
"""
pow5bits(e) = ((e * 1217359) >> 19) + 1

""""
Ryu.mulshift(m::UInt64, mula, mulb, j)::UInt64
Compute `(m * mul) >> j`, where `mul = mula + mulb<<64` and `j >= 64`.
"""
@inline mulshift(m::UInt64, mula, mulb, j) = ((((UInt128(m) * mula) >> 64) + UInt128(m) * mulb) >> (j - 64)) % UInt64

""""
Ryu.mulshift(m::UInt32, mul, j)::UInt32
Compute `(m * mul) >> j`, where `j >= 32`.
"""
@inline mulshift(m::UInt32, mul, j) = ((((UInt64(m) * (mul % UInt32)) >> 32) + (UInt64(m) * (mul >> 32))) >> (j - 32)) % UInt32

""""
Ryu.mulshift(m::UInt16, mul, j)::UInt16
Compute `(m * mul) >> j`, where `j >= 16`.
"""
@inline mulshift(m::UInt16, mul, j) = ((((UInt32(m) * (mul % UInt16)) >> 16) + (UInt32(m) * (mul >> 16))) >> (j - 16))

indexforexp(e) = div(e + 15, 16)
pow10bitsforindex(idx) = 16 * idx + 120
lengthforindex(idx) = div(((Int64(16 * idx) * 1292913986) >> 32) + 1 + 16 + 8, 9)

"""
Ryu.pow5(x, p)
Return `true` if `5^p` is a divisor of `x`.
"""
@inline function pow5(x, p)
count = 0
while true
Expand All @@ -57,8 +83,18 @@ lengthforindex(idx) = div(((Int64(16 * idx) * 1292913986) >> 32) + 1 + 16 + 8, 9
end
end

"""
Ryu.pow2(x, p)
Return `true` if `2^p` is a divisor of `x`. In other words, if the trailing `p` bits of `x` are zero.
"""
pow2(x, p) = (x & ((Int64(1) << p) - 1)) == 0

"""
Ryu.decimallength(v)
The number of decimal digits of the integer `v`.
"""
@inline function decimallength(v)
v >= 10000000000000000 && return 17
v >= 1000000000000000 && return 16
Expand All @@ -78,7 +114,6 @@ pow2(x, p) = (x & ((Int64(1) << p) - 1)) == 0
v >= 10 && return 2
return 1
end

@inline function decimallength(v::UInt32)
v >= 100000000 && return 9
v >= 10000000 && return 8
Expand All @@ -90,7 +125,6 @@ end
v >= 10 && return 2
return 1
end

@inline function decimallength(v::UInt16)
v >= 10000 && return 5
v >= 1000 && return 4
Expand Down Expand Up @@ -147,6 +181,11 @@ end
return vr, vp, vm
end

"""
Ryu.umul256(a::UInt128, bHi::UInt64, bLo::UInt64)::Tuple{UInt128, UInt128}
Compute `p = a*b` where `b = bLo + bHi<<64`, returning the result as `pLo, pHi` where `p = pLo + pHi<<128`.
"""
@inline function umul256(a, bHi, bLo)
aLo = a % UInt64
aHi = (a >> 64) % UInt64
Expand All @@ -172,8 +211,18 @@ end
return pLo, pHi
end

"""
Ryu.umul256_hi(a::UInt128, bHi::UInt64, bLo::UInt64)::UInt128
Compute `pHi = (a*b)>>128` where `b = bLo + bHi<<64`.
"""
@inline umul256_hi(a, bHi, bLo) = umul256(a, bHi, bLo)[2]

"""
Ryu.mulshiftmod1e9(m, mula, mulb, mulc, j)::UInt32
Compute `(m * mul) >> j % 10^9` where `mul = mula + mulb<<64 + mulc<<128`, and `j >= 128`.
"""
@inline function mulshiftmod1e9(m, mula, mulb, mulc, j)
b0 = UInt128(m) * mula
b1 = UInt128(m) * mulb
Expand Down Expand Up @@ -323,40 +372,38 @@ end

const POW10_OFFSET_2, MIN_BLOCK_2, POW10_SPLIT_2 = generateinversetables()

bitlength(this) = Base.GMP.MPZ.sizeinbase(this, 2)

@inline function pow5invsplit(::Type{Float64}, i)
pow = big(5)^i
inv = div(big(1) << (bitlength(pow) - 1 + pow5_inv_bitcount(Float64)), pow) + 1
inv = div(big(1) << (ndigits(pow, base=2) - 1 + pow5_inv_bitcount(Float64)), pow) + 1
return (UInt64(inv & ((big(1) << 64) - 1)), UInt64(inv >> 64))
end

@inline function pow5invsplit(::Type{Float32}, i)
pow = big(5)^i
inv = div(big(1) << (bitlength(pow) - 1 + pow5_inv_bitcount(Float32)), pow) + 1
inv = div(big(1) << (ndigits(pow, base=2) - 1 + pow5_inv_bitcount(Float32)), pow) + 1
return UInt64(inv)
end

@inline function pow5invsplit(::Type{Float16}, i)
pow = big(5)^i
inv = div(big(1) << (bitlength(pow) - 1 + pow5_inv_bitcount(Float16)), pow) + 1
inv = div(big(1) << (ndigits(pow, base=2) - 1 + pow5_inv_bitcount(Float16)), pow) + 1
return UInt32(inv)
end

@inline function pow5split(::Type{Float64}, i)
pow = big(5)^i
j = bitlength(pow) - pow5_bitcount(Float64)
j = ndigits(pow, base=2) - pow5_bitcount(Float64)
return (UInt64((pow >> j) & ((big(1) << 64) - 1)), UInt64(pow >> (j + 64)))
end

@inline function pow5split(::Type{Float32}, i)
pow = big(5)^i
return UInt64(pow >> (bitlength(pow) - pow5_bitcount(Float32)))
return UInt64(pow >> (ndigits(pow, base=2) - pow5_bitcount(Float32)))
end

@inline function pow5split(::Type{Float16}, i)
pow = big(5)^i
return UInt32(pow >> (bitlength(pow) - pow5_bitcount(Float16)))
return UInt32(pow >> (ndigits(pow, base=2) - pow5_bitcount(Float16)))
end

const DOUBLE_POW5_INV_SPLIT = map(i->pow5invsplit(Float64, i), 0:291)
Expand Down

0 comments on commit d036575

Please sign in to comment.