Skip to content

Commit

Permalink
Bug fixes to the vectorized versions of a few hyp and trig functions.
Browse files Browse the repository at this point in the history
  • Loading branch information
chriselrod committed Dec 25, 2018
1 parent 0cb1736 commit 8b83a5a
Show file tree
Hide file tree
Showing 3 changed files with 30 additions and 25 deletions.
22 changes: 11 additions & 11 deletions src/hyp.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,12 +8,12 @@ over_sch(::Type{Float32}) = 89f0
Compute hyperbolic sine of `x`.
"""
function sinh(x::FloatType)
function sinh(x::V) where V <: FloatType
T = eltype(x)
u = abs(x)
d = expk2(Double(u))
d = dsub(d, drec(d))
u = T(d) * T(0.5)
u = V(d) * T(0.5)
u = vifelse(abs(x) > over_sch(T), T(Inf), u)
u = vifelse(isnan(u), T(Inf), u)
u = flipsign(u, x)
Expand All @@ -28,12 +28,12 @@ end
Compute hyperbolic cosine of `x`.
"""
function cosh(x::FloatType)
function cosh(x::V) where V <: FloatType
T = eltype(x)
u = abs(x)
d = expk2(Double(u))
d = dadd(d, drec(d))
u = T(d) * T(0.5)
u = V(d) * T(0.5)
u = vifelse(abs(x) > over_sch(T), T(Inf), u)
u = vifelse(isnan(u), T(Inf), u)
u = vifelse(isnan(x), T(NaN), u)
Expand All @@ -50,13 +50,13 @@ over_th(::Type{Float32}) = 18.714973875f0
Compute hyperbolic tangent of `x`.
"""
function tanh(x::FloatType)
function tanh(x::V) where V <: FloatType
T = eltype(x)
u = abs(x)
d = expk2(Double(u))
e = drec(d)
d = ddiv(dsub(d, e), dadd(d, e))
u = T(d)
u = V(d)
u = vifelse(abs(x) > over_th(T), T(1.0), u)
u = vifelse(isnan(u), T(1), u)
u = flipsign(u, x)
Expand All @@ -81,7 +81,7 @@ function asinh(x::V) where V <: FloatType
d = vifelse(yg1, dmul(d, y), d)

d = logk2(dnormalize(dadd(d, x)))
y = T(d)
y = V(d)

y = vifelse(((abs(x) > SQRT_MAX(T)) | isnan(y)), flipsign(T(Inf), x), y)
y = vifelse(isnan(x), T(NaN), y)
Expand All @@ -97,10 +97,10 @@ end
Compute the inverse hyperbolic cosine of `x`.
"""
function acosh(x::FloatType)
function acosh(x::V) where V <: FloatType
T = eltype(x)
d = logk2(dadd2(dmul(dsqrt(dadd2(x, T(1.0))), dsqrt(dsub2(x, T(1.0)))), x))
y = T(d)
y = V(d)

y = vifelse(((x > SQRT_MAX(T)) | isnan(y)), T(Inf), y)
y = vifelse(x == T(1.0), T(0.0), y)
Expand All @@ -117,11 +117,11 @@ end
Compute the inverse hyperbolic tangent of `x`.
"""
function atanh(x::FloatType)
function atanh(x::V) where V <: FloatType
T = eltype(x)
u = abs(x)
d = logk2(ddiv(dadd2(T(1.0), u), dsub2(T(1.0), u)))
u = vifelse(u > T(1.0), T(NaN), vifelse(u == T(1.0), T(Inf), T(d) * T(0.5)))
u = vifelse(u > T(1.0), T(NaN), vifelse(u == T(1.0), T(Inf), V(d) * T(0.5)))

u = vifelse(isinf(x) | isnan(u), T(NaN), u)
u = flipsign(u, x)
Expand Down
14 changes: 9 additions & 5 deletions src/priv.jl
Original file line number Diff line number Diff line change
Expand Up @@ -150,7 +150,11 @@ end

@inline function atan2k(y::Double{V}, x::Double{V}) where {T,V<:Union{T,Vec{<:Any,T}}}
xl0 = x < 0
q = vifelse(xl0, -2, 0)
if V <: Vec
q = vifelse(xl0, Vec{length(x.hi),EquivalentInteger(T)}(-2), 0)
else
q = vifelse(xl0, -2, 0)
end
x = vifelse(xl0, -x, x)
ygx = y > x
t = x
Expand All @@ -167,7 +171,7 @@ end
t = dmul(t, u)
t = dmul(s, dadd(T(1.0), t))
T <: Float64 && (t = vifelse(abs(s.hi) < 1e-200, s, t))
t = dadd(dmul(V(q), MDPI2(T)), t)
t = dadd(dmul(convert(V,q), MDPI2(T)), t)
return t
end

Expand Down Expand Up @@ -247,8 +251,8 @@ end
end

@inline function expk2(d::Double{V}) where {T<:Union{Float32,Float64},V<:Union{T,Vec{<:Any,T}}}
q = round(T(d) * T(MLN2E))
qi = unsafe_trunc(Int, q)
q = round(V(d) * T(MLN2E))
qi = unsafe_trunc(EquivalentInteger(T), q)

s = dadd(d, -q * L2U(T))
s = dadd(s, -q * L2L(T))
Expand Down Expand Up @@ -295,7 +299,7 @@ end

t = logk2_kernel(x2.hi)

s = dmul(MDLN2(T), T(e))
s = dmul(MDLN2(T), convert(V,e))
s = dadd(s, scale(x, T(2.0)))
s = dadd(s, dmul(dmul(x2, x), t))
return s
Expand Down
19 changes: 10 additions & 9 deletions src/trig.jl
Original file line number Diff line number Diff line change
Expand Up @@ -650,7 +650,7 @@ function tan_fast(d::FloatType32)
end


@inline function tan_kernel(x::Double{Float64})
@inline function tan_kernel(x::Double{<:FloatType64})
c15 = 1.01419718511083373224408e-05
c14 = -2.59519791585924697698614e-05
c13 = 5.23388081915899855325186e-05
Expand Down Expand Up @@ -680,11 +680,11 @@ end
return dadd(c1, x.hi * (@horner x.hi c2 c3 c4 c5 c6 c7))
end

function tan(d::FloatType64)
function tan(d::V) where V <: FloatType64
T = eltype(d)
qh = trunc(d * (T(M_2_PI)) / (1 << 24))
s = dadd2(dmul(Double(T(M_2_PI_H), T(M_2_PI_L)), d), (d < 0 ? T(-0.5) : T(0.5)) - qh * (1 << 24))
ql = trunc(T(s))
s = dadd2(dmul(Double(T(M_2_PI_H), T(M_2_PI_L)), d), vifelse(d < 0, V(-0.5), V(0.5)) - qh * (1 << 24))
ql = trunc(V(s))

s = dadd2(d, qh * (-PI_A(T) * T(0.5) * (1 << 24)))
s = dadd2(s, ql * (-PI_A(T) * T(0.5) ))
Expand All @@ -708,7 +708,7 @@ function tan(d::FloatType64)

x = vifelse(qli_odd, drec(x), x)

v = T(x)
v = V(x)

v = vifelse(
(!isinf(d)) & (isnegzero(d) | (abs(d) > TRIG_MAX(T))),
Expand All @@ -718,7 +718,7 @@ function tan(d::FloatType64)
return v
end

function tan(d::FloatType32)
function tan(d::V) where V <: FloatType32
T = eltype(d)
q = round(d * (T(M_2_PI)))

Expand All @@ -743,7 +743,7 @@ function tan(d::FloatType32)

x = vifelse(qli_odd, drec(x), x)

v = T(x)
v = V(x)

v = vifelse(
(!isinf(d)) & (isnegzero(d) | (abs(d) > TRIG_MAX(T))),
Expand Down Expand Up @@ -900,12 +900,13 @@ end
Compute the inverse cosine of `x`, where the output is in radians.
"""
function acos(x::T) where {T<:FloatType}
function acos(x::V) where {V<:FloatType}
T = eltype(x)
d = atan2k(dsqrt(dmul(dadd(T(1), x), dsub(T(1), x))), Double(abs(x)))
d = flipsign(d, x)
d = vifelse(abs(x) == 1, Double(T(0)), d)
d = vifelse(signbit(x), dadd(MDPI(T), d), d)
return T(d)
return V(d)
end


Expand Down

0 comments on commit 8b83a5a

Please sign in to comment.