Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fast exp2(::Int) (fixs #17412) #17447

Merged
merged 10 commits into from
Aug 1, 2016
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions base/float.jl
Original file line number Diff line number Diff line change
Expand Up @@ -519,9 +519,9 @@ exponent_one(::Type{Float32}) = 0x3f80_0000
exponent_half(::Type{Float32}) = 0x3f00_0000
significand_mask(::Type{Float32}) = 0x007f_ffff

significand_bits{T<:AbstractFloat}(::Type{T}) = trailing_ones(significand_mask(T))
exponent_bits{T<:AbstractFloat}(::Type{T}) = sizeof(T)*8 - significand_bits(T) - 1
exponent_bias{T<:AbstractFloat}(::Type{T}) = Int(exponent_one(T) >> significand_bits(T))
@pure significand_bits{T<:AbstractFloat}(::Type{T}) = trailing_ones(significand_mask(T))
@pure exponent_bits{T<:AbstractFloat}(::Type{T}) = sizeof(T)*8 - significand_bits(T) - 1
@pure exponent_bias{T<:AbstractFloat}(::Type{T}) = Int(exponent_one(T) >> significand_bits(T))

## Array operations on floating point numbers ##

Expand Down
15 changes: 15 additions & 0 deletions base/math.jl
Original file line number Diff line number Diff line change
Expand Up @@ -124,6 +124,21 @@ for f in (:sinh, :cosh, :tanh, :atan, :asinh, :exp, :erf, :erfc, :expm1)
@eval ($f)(x::AbstractFloat) = error("not implemented for ", typeof(x))
end

# functions with special cases for integer arguments
@inline function exp2(x::Base.BitInteger)
if x > 1023
Inf64
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

should exp2(x::BigInt) return a BigFloat answer?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Nice catch.
I've added some test cases to test/bigint.jl to ensure this.

elseif x <= -1023
# if -1073 < x <= -1023 then Result will be a subnormal number
# Hex literal with padding must be use to work on 32bit machine
reinterpret(Float64, 0x0000_0000_0000_0001 << ((x + 1074)) % UInt)
else
# We will cast everything to Int64 to avoid errors in case of Int128
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Since you cast everything to Int64, you could also restrict the argument type to Int, and have a generic version that accepts Integer and calls the main function with argument converted to Int.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think it is a bit cleaner not to.
Not all paths have the cast, on x
only line 141

The cast on the 1 in line 134, is to make sure things work on 32 bit machines where 1 is a 32bit.
Perhaps cleaner, I change that one to 0x0000_0000_0000_0001

# If x is a Int128, and is outside the range of Int64, then it is not -1023<x<=1023
reinterpret(Float64, (exponent_bias(Float64) + (x % Int64)) << (significand_bits(Float64)) % UInt)
end
end

# TODO: GNU libc has exp10 as an extension; should openlibm?
exp10(x::Float64) = 10.0^x
exp10(x::Float32) = 10.0f0^x
Expand Down
19 changes: 19 additions & 0 deletions test/bigint.jl
Original file line number Diff line number Diff line change
Expand Up @@ -311,3 +311,22 @@ ndigits(rand(big(-999:999)), big(2)^rand(2:999))

@test big(5)^true == big(5)
@test big(5)^false == one(BigInt)


# operations that when applied to Int64 give Float64, should give BigFloat
@test typeof(exp(a)) == BigFloat
@test typeof(exp2(a)) == BigFloat
@test typeof(exp10(a)) == BigFloat
@test typeof(expm1(a)) == BigFloat
@test typeof(erf(a)) == BigFloat
@test typeof(erfc(a)) == BigFloat
@test typeof(cosh(a)) == BigFloat
@test typeof(sinh(a)) == BigFloat
@test typeof(tanh(a)) == BigFloat
@test typeof(sech(a)) == BigFloat
@test typeof(csch(a)) == BigFloat
@test typeof(coth(a)) == BigFloat
@test typeof(cbrt(a)) == BigFloat
@test typeof(tan(a)) == BigFloat
@test typeof(cos(a)) == BigFloat
@test typeof(sin(a)) == BigFloat
16 changes: 16 additions & 0 deletions test/math.jl
Original file line number Diff line number Diff line change
Expand Up @@ -192,6 +192,22 @@ end
@test exp2(Float16(2.)) ≈ exp2(2.)
@test log(e) == 1

# check exp2(::Integer) matches exp2(::Float)"
for ii in -2048:2048
expected = exp2(float(ii))
@test(exp2(Int16(ii)) == expected)
@test(exp2(Int32(ii)) == expected)
@test(exp2(Int64(ii)) == expected)
@test(exp2(Int128(ii)) == expected)
if ii >= 0
@test(exp2(UInt16(ii)) == expected)
@test(exp2(UInt32(ii)) == expected)
@test(exp2(UInt64(ii)) == expected)
@test(exp2(UInt128(ii)) == expected)
end
end


for T in (Int, Float64, BigFloat)
@test deg2rad(T(180)) ≈ 1pi
@test deg2rad(T[45, 60]) ≈ [pi/T(4), pi/T(3)]
Expand Down