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

Add Julia ports of hyperbolic (including inverse) functions from openlibm. #24159

Merged
merged 28 commits into from
Feb 14, 2018
Merged
Show file tree
Hide file tree
Changes from 1 commit
Commits
Show all changes
28 commits
Select commit Hold shift + click to select a range
fd37f3d
Add Julia ports of hyperbolic functions from openlibm.
pkofod Oct 13, 2017
0bead47
Actually save the correct tests...
pkofod Oct 16, 2017
53886c5
Fix comments in hyperbolic functions and add acosh.
pkofod Oct 16, 2017
9cba2fa
Use Julia version of ldexp instead of openlibm.
pkofod Oct 16, 2017
423d46a
Add asinh
pkofod Oct 17, 2017
dd8c892
Fix missing end and improve comments.
pkofod Oct 17, 2017
4daaa8b
Spelling.
pkofod Oct 17, 2017
6bd838c
Changes
pkofod Oct 17, 2017
5f8f975
Add atanh and minor fixes in general.
pkofod Oct 17, 2017
ecbc401
Remove atanh from loop in math.jl.
pkofod Oct 17, 2017
2aec485
_ldexp_exp to ldexp
pkofod Oct 18, 2017
7bf4cdf
Readd ldexp_exp
pkofod Nov 23, 2017
90b5c45
Remove semi-colon.
pkofod Nov 23, 2017
ea5c860
Add tests for inverse hyperbolics.
pkofod Nov 23, 2017
79d7477
Fix licences.
pkofod Nov 23, 2017
ee9b458
Fix a missing f0.
pkofod Nov 24, 2017
1fbca61
Clean up tested values.
pkofod Nov 24, 2017
0424b85
Merge branch 'master' into hyperbolic
simonbyrne Jan 28, 2018
cc7e1f7
Fix atanh for NaN input.
pkofod Feb 6, 2018
248dcbb
Whitespace.
pkofod Feb 6, 2018
2b396ad
Merge branch 'master' into hyperbolic
pkofod Feb 12, 2018
c15347b
Improve domain error messages.
pkofod Feb 12, 2018
79ae7d1
Make tests neater.
pkofod Feb 12, 2018
eb1390a
Merge branch 'hyperbolic' of github.com:pkofod/julia into hyperbolic
pkofod Feb 12, 2018
c3bedbd
Fix wrong merge conflict resolution, and make tests prettier.
pkofod Feb 13, 2018
d92b7bc
More testing.
pkofod Feb 13, 2018
b751c37
Fix end
pkofod Feb 13, 2018
5c22157
More test changes.
pkofod Feb 13, 2018
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
2 changes: 1 addition & 1 deletion base/math.jl
Original file line number Diff line number Diff line change
Expand Up @@ -244,7 +244,7 @@ asinh(x::Number)
Accurately compute ``e^x-1``.
"""
expm1(x)
for f in (:cbrt, :atan, :exp2, :expm1)
for f in (:cbrt, :exp2, :expm1)
@eval begin
($f)(x::Float64) = ccall(($(string(f)),libm), Float64, (Float64,), x)
($f)(x::Float32) = ccall(($(string(f,"f")),libm), Float32, (Float32,), x)
Expand Down
18 changes: 8 additions & 10 deletions base/special/hyperbolic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -18,10 +18,10 @@ _ldexp_exp(x::Real, i) = _ldexp_exp(float(x, i))

# Hyperbolic functions
# sinh methods
SINH_SMALL_X(::Type{Float64}) = 2.0^-28
H_SMALL_X(::Type{Float64}) = 2.0^-28
H_MEDIUM_X(::Type{Float64}) = 22.0

SINH_SMALL_X(::Type{Float32}) = 2f-12
H_SMALL_X(::Type{Float32}) = 2f-12
H_MEDIUM_X(::Type{Float32}) = 9f0

H_LARGE_X(::Type{Float64}) = 709.7822265633563 # nextfloat(709.7822265633562)
Expand All @@ -34,9 +34,9 @@ function sinh(x::T) where T <: Union{Float32, Float64}
# mathematically sinh(x) is defined to be (exp(x)-exp(-x))/2
# 1. Replace x by |x| (sinh(-x) = -sinh(x)).
# 2. Find the branch and the expression to calculate and return it
# a) 0 <= x < SINH_SMALL_X
# a) 0 <= x < H_SMALL_X
# return x
# b) SINH_SMALL_X <= x < H_MEDIUM_X
# b) H_SMALL_X <= x < H_MEDIUM_X
# return sinh(x) = (E + E/(E+1))/2, where E=expm1(x)
# c) H_MEDIUM_X <= x < H_LARGE_X
# return sinh(x) = exp(x)/2
Expand All @@ -59,7 +59,7 @@ function sinh(x::T) where T <: Union{Float32, Float64}
# in a) or b)
if absx < H_MEDIUM_X(T)
# in a)
if absx < SINH_SMALL_X(T)
if absx < H_SMALL_X(T)
return x
end
t = expm1(absx)
Expand Down Expand Up @@ -137,16 +137,14 @@ cosh(x::Real) = cosh(float(x))
# tanh methods
TANH_LARGE_X(::Type{Float64}) = 22.0
TANH_LARGE_X(::Type{Float32}) = 9.0f0
TANH_SMALL_X(::Type{Float64}) = 2.0^-28
TANH_SMALL_X(::Type{Float32}) = 2f0^-12
function tanh(x::T) where T<:Union{Float32, Float64}
# Method
# mathematically tanh(x) is defined to be (exp(x)-exp(-x))/(exp(x)+exp(-x))
# 1. reduce x to non-negative by tanh(-x) = -tanh(x).
# 2. Find the branch and the expression to calculate and return it
# a) 0 <= x < TANH_SMALL_X
# a) 0 <= x < H_SMALL_X
# return x
# b) TANH_SMALL_X <= x < 1
# b) H_SMALL_X <= x < 1
# -expm1(-2x)/(expm1(-2x) + 2)
# c) 1 <= x < TANH_LARGE_X
# 1 - 2/(expm1(2x) + 2)
Expand All @@ -161,7 +159,7 @@ function tanh(x::T) where T<:Union{Float32, Float64}
absx = abs(x)
if absx < TANH_LARGE_X(T)
# in a)
if absx < TANH_SMALL_X(T)
if absx < H_SMALL_X(T)
return x
end
if absx >= T(1)
Expand Down
34 changes: 14 additions & 20 deletions test/math.jl
Original file line number Diff line number Diff line change
Expand Up @@ -844,6 +844,8 @@ end
end
#prev, current, next float
pcnfloat(x) = prevfloat(x), x, nextfloat(x)
import Base.Math: COSH_SMALL_X, H_SMALL_X, H_MEDIUM_X, H_LARGE_X

@testset "sinh" begin
for T in (Float32, Float64)
@test sinh(zero(T)) === zero(T)
Expand All @@ -853,14 +855,10 @@ pcnfloat(x) = prevfloat(x), x, nextfloat(x)
@test sinh(T(1000)) === T(Inf)
@test sinh(-T(1000)) === -T(Inf)
@test isnan_type(T, sinh(T(NaN)))
end
for x in (pcnfloat(2.0^-28)..., pcnfloat(22.0)..., pcnfloat(709.7822265633563)...)
@test sinh(x) ≈ sinh(big(x)) rtol=eps(Float64)
@test sinh(-x) ≈ sinh(big(-x)) rtol=eps(Float64)
end
for x in (pcnfloat(2f-12)..., pcnfloat(9f0)..., pcnfloat(88.72283f0)...)
@test sinh(x) ≈ sinh(big(x)) rtol=eps(Float32)
@test sinh(-x) ≈ sinh(big(-x)) rtol=eps(Float32)
for x in Iterators.flatten(pcnfloat.([H_SMALL_X(T), H_MEDIUM_X(T), H_LARGE_X(T)]))
@test sinh(x) ≈ sinh(big(x)) rtol=eps(T)
@test sinh(-x) ≈ sinh(big(-x)) rtol=eps(T)
end
end
end

Expand All @@ -873,14 +871,10 @@ end
@test cosh(T(1000)) === T(Inf)
@test cosh(-T(1000)) === T(Inf)
@test isnan_type(T, cosh(T(NaN)))
end
for x in (pcnfloat(2.7755602085408512e-17)..., pcnfloat(22.0)..., pcnfloat(709.7822265633563)...)
@test cosh(x) ≈ cosh(big(x)) rtol=eps(Float64)
@test cosh(-x) ≈ cosh(big(-x)) rtol=eps(Float64)
end
for x in (pcnfloat(0.00024414062f0)..., pcnfloat(9f0)..., pcnfloat(88.72283f0)...)
@test cosh(x) ≈ cosh(big(x)) rtol=eps(Float32)
@test cosh(-x) ≈ cosh(big(-x)) rtol=eps(Float32)
for x in Iterators.flatten(pcnfloat.([COSH_SMALL_X(T), H_MEDIUM_X(T), H_LARGE_X(T)]))
@test cosh(x) ≈ cosh(big(x)) rtol=eps(T)
@test cosh(-x) ≈ cosh(big(-x)) rtol=eps(T)
end
end
end

Expand All @@ -893,7 +887,7 @@ end
@test tanh(T(1000)) === one(T)
@test tanh(-T(1000)) === -one(T)
@test isnan_type(T, tanh(T(NaN)))
for x in Iterators.flatten(pcnfloat.(T == Float64 ? [2.0^-28,1.0,22.0] : [2f0^-12,1f0,9f0]))
for x in Iterators.flatten(pcnfloat.([H_SMALL_X(T), T(1.0), H_MEDIUM_X(T)]))
@test tanh(x) ≈ tanh(big(x)) rtol=eps(T)
@test tanh(-x) ≈ tanh(big(-x)) rtol=eps(T)
end
Expand All @@ -907,7 +901,7 @@ end
@test asinh(nextfloat(zero(T))) === nextfloat(zero(T))
@test asinh(prevfloat(zero(T))) === prevfloat(zero(T))
@test isnan_type(T, asinh(T(NaN)))
for x in Iterators.flatten(pcnfloat.(T == Float64 ? [2.0^-28,2.0,2.0^28] : [2f0^-12,2f0,22f28]))
for x in Iterators.flatten(pcnfloat.([T(2)^-28,T(2),T(2)]))
Copy link
Member

Choose a reason for hiding this comment

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

The second T(2) doesn't look right. I assume you mean T(2)^28?

@test asinh(x) ≈ asinh(big(x)) rtol=eps(T)
@test asinh(-x) ≈ asinh(big(-x)) rtol=eps(T)
end
Expand All @@ -919,7 +913,7 @@ end
@test_throws DomainError acosh(T(0.1))
@test acosh(one(T)) === zero(T)
@test isnan_type(T, acosh(T(NaN)))
for x in (nextfloat(T(1)), pcnfloat(T(2))..., pcnfloat(T(2^28))...)
for x in Iterators.flatten(pcnfloat.([T(1) T(2), T(2^28)]))
Copy link
Member

@stevengj stevengj Feb 13, 2018

Choose a reason for hiding this comment

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

Shouldn't it be Iterators.flatten([nextfloat(T(1)); pcnfloat.([T(2), T(2^28)])])? Or [nextfloat(T(1)), pcnfloat(T(2))..., pcnfloat(T(2^28))...]

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 seem to have pushed without this fix.

@test acosh(x) ≈ acosh(big(x) rtol=eps(T)
end
end
Expand All @@ -935,7 +929,7 @@ end
@test atanh(nextfloat(zero(T))) === nextfloat(zero(T))
@test atanh(prevfloat(zero(T))) === prevfloat(zero(T))
@test isnan_type(T, atanh(T(NaN)))
for x in Iterators.flatten(pcnfloat.(T == Float64 ? [2.0^-28,0.5] : [2f-28,0.5f0]))
for x in Iterators.flatten(pcnfloat.([T(2.0)^-28, T(0.5)]))
@test atanh(x) ≈ atanh(big(x)) rtol=eps(T)
@test atanh(-x) ≈ atanh(big(-x)) rtol=eps(T)
end
Expand Down