Merge pull request #204 from araujoms/main

Fix inf and nan bugs
JuliaMath · Jun 27, 2024 · 47c49f0 · 47c49f0
2 parents 9e9d5b3 + 60e845c
commit 47c49f0
Show file tree

Hide file tree

Showing 9 changed files with 50 additions and 34 deletions.
diff --git a/src/math/elementary/explog.jl b/src/math/elementary/explog.jl
@@ -259,7 +259,7 @@ end
 
 function log1p(x::DoubleFloat{T})  where {T<:IEEEFloat}
     isnan(x) && return x
-    isinf(x) && !signbit(x) && return 
+    isinf(x) && !signbit(x) && return x
     u = 1.0 + x
     if u == one(DoubleFloat{T})
         x

diff --git a/src/math/ops.jl b/src/math/ops.jl
@@ -6,7 +6,7 @@ isnotfinite(x::Float32) = reinterpret(UInt32,x) & 0x7f800000 === 0x7f800000
 isqnan(x::Float32) = reinterpret(UInt32,x)      & 0x7fc00000 === 0x7fc00000
 isainf(x::Float32) = reinterpret(UInt32,x)      & 0x7fc00000 === 0x7f800000
 
-isnotfinite(x::DoubleFloat{T}) where {T<:IEEEFloat} = isnan(LO(x))
+isnotfinite(x::DoubleFloat{T}) where {T<:IEEEFloat} = !isfinite(HI(x))
 isqnan(x::DoubleFloat{T}) where {T<:IEEEFloat} = isqnan(HI(x))
 isainf(x::DoubleFloat{T}) where {T<:IEEEFloat} = isainf(HI(x))
 

diff --git a/src/math/ops/op_dbdb_db.jl b/src/math/ops/op_dbdb_db.jl
@@ -12,41 +12,41 @@ Base.:(/)(x::Tuple{T,T}, y::DoubleFloat{T}) where {T<:IEEEFloat} = DoubleFloat{T
 
 
 @inline function add_dbdb_db(x::DoubleFloat{T}, y::DoubleFloat{T}) where {T<:IEEEFloat}
-    (isnan(LO(x)) | isnan(LO(y))) && return add_dbdb_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return add_dbdb_db_nonfinite(x,y)
     return DoubleFloat{T}(add_dddd_dd(HILO(x), HILO(y)))
 end
 
 @inline function sub_dbdb_db(x::DoubleFloat{T}, y::DoubleFloat{T}) where {T<:IEEEFloat}
-    (isnan(LO(x)) | isnan(LO(y))) && return sub_dbdb_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return sub_dbdb_db_nonfinite(x,y)
     return DoubleFloat{T}(sub_dddd_dd(HILO(x), HILO(y)))
 end
 
 @inline function mul_dbdb_db(x::DoubleFloat{T}, y::DoubleFloat{T}) where {T<:IEEEFloat}
-    (isnan(LO(x)) | isnan(LO(y))) && return mul_dbdb_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return mul_dbdb_db_nonfinite(x,y)
     return DoubleFloat{T}(mul_dddd_dd(HILO(x), HILO(y)))
 end
 
 @inline function dvi_dbdb_db(x::DoubleFloat{T}, y::DoubleFloat{T}) where {T<:IEEEFloat}
-    (iszero(HI(y)) | isnan(LO(x)) | isnan(LO(y))) && return dvi_dbdb_db_nonfinite(x,y)
+    (iszero(HI(y)) | isnotfinite(x) | isnotfinite(y)) && return dvi_dbdb_db_nonfinite(x,y)
     return DoubleFloat{T}(dvi_dddd_dd(HILO(x), HILO(y)))
 end
 
 @inline function add_dbdb_db_nonfinite(x::DoubleFloat{T}, y::DoubleFloat{T}) where {T<:IEEEFloat}
     z = HI(x) + HI(y)
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function sub_dbdb_db_nonfinite(x::DoubleFloat{T}, y::DoubleFloat{T}) where {T<:IEEEFloat}
     z = HI(x) - HI(y)
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function mul_dbdb_db_nonfinite(x::DoubleFloat{T}, y::DoubleFloat{T}) where {T<:IEEEFloat}
     z = HI(x) * HI(y)
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function dvi_dbdb_db_nonfinite(x::DoubleFloat{T}, y::DoubleFloat{T}) where {T<:IEEEFloat}
     z = HI(x) / HI(y)
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
diff --git a/src/math/ops/op_dbfp_db.jl b/src/math/ops/op_dbfp_db.jl
@@ -1,39 +1,39 @@
 @inline function add_dbfp_db(x::DoubleFloat{T}, y::T) where {T<:IEEEFloat}
-    (isqnan(LO(x)) | isnotfinite(y)) && return add_dbfp_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return add_dbfp_db_nonfinite(x,y)
     return DoubleFloat{T}(add_ddfp_dd(HILO(x), y))
 end
 
 @inline function sub_dbfp_db(x::DoubleFloat{T}, y::T) where {T<:IEEEFloat}
-    (isqnan(LO(x)) | isnotfinite(y)) && return sub_dbfp_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return sub_dbfp_db_nonfinite(x,y)
     return DoubleFloat{T}(sub_ddfp_dd(HILO(x), y))
 end
 
 @inline function mul_dbfp_db(x::DoubleFloat{T}, y::T) where {T<:IEEEFloat}
-    (isqnan(LO(x)) | isnotfinite(y)) && return mul_dbfp_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return mul_dbfp_db_nonfinite(x,y)
     return DoubleFloat{T}(mul_ddfp_dd(HILO(x), y))
 end
 
 @inline function dvi_dbfp_db(x::DoubleFloat{T}, y::T) where {T<:IEEEFloat}
-    (isqnan(LO(x)) | isnotfinite(y)) && return dvi_dbfp_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return dvi_dbfp_db_nonfinite(x,y)
     return DoubleFloat{T}(dvi_ddfp_dd(HILO(x), y))
 end
 
 @inline function add_dbfp_db_nonfinite(x::DoubleFloat{T}, y::T) where {T<:IEEEFloat}
     z = HI(x) + y
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function sub_dbfp_db_nonfinite(x::DoubleFloat{T}, y::T) where {T<:IEEEFloat}
     z = HI(x) - y
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function mul_dbfp_db_nonfinite(x::DoubleFloat{T}, y::T) where {T<:IEEEFloat}
     z = HI(x) * y
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function dvi_dbfp_db_nonfinite(x::DoubleFloat{T}, y::T) where {T<:IEEEFloat}
     z = HI(x) / y
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
diff --git a/src/math/ops/op_dd_dd.jl b/src/math/ops/op_dd_dd.jl
@@ -56,8 +56,9 @@ end
 end
 
 function sqrt_dd_dd(x::Tuple{T,T}) where {T<:IEEEFloat}
-    iszero(HI(x)) && return x
-    signbit(HI(x)) && throw(DomainError("sqrt(x) expects x >= 0"))
+    (isnan(HI(x)) | iszero(HI(x))) && return x
+    signbit(HI(x)) && Base.Math.throw_complex_domainerror(:sqrt, HI(x))
+    isinf(HI(x)) && return x
 
     ahi, alo = HILO(x)
     s = sqrt(ahi)

diff --git a/src/math/ops/op_fpdb_db.jl b/src/math/ops/op_fpdb_db.jl
@@ -1,39 +1,39 @@
 @inline function add_fpdb_db(x::T, y::DoubleFloat{T}) where {T<:IEEEFloat}
-    (isnotfinite(x) | isqnan(LO(y))) && return add_fpdb_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return add_fpdb_db_nonfinite(x,y)
     return DoubleFloat{T}(add_fpdd_dd(x, HILO(y)))
 end
 
 @inline function sub_fpdb_db(x::T, y::DoubleFloat{T}) where {T<:IEEEFloat}
-    (isnotfinite(x) | isqnan(LO(y))) && return sub_fpdb_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return sub_fpdb_db_nonfinite(x,y)
     return DoubleFloat{T}(sub_fpdd_dd(x, HILO(y)))
 end
 
 @inline function mul_fpdb_db(x::T, y::DoubleFloat{T}) where {T<:IEEEFloat}
-    (isnotfinite(x) | isqnan(LO(y))) && return mul_fpdb_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return mul_fpdb_db_nonfinite(x,y)
     return DoubleFloat{T}(mul_fpdd_dd(x, HILO(y)))
 end
 
 @inline function dvi_fpdb_db(x::T, y::DoubleFloat{T}) where {T<:IEEEFloat}
-    (isnotfinite(x) | isqnan(LO(y))) && return dvi_fpdb_db_nonfinite(x,y)
+    (isnotfinite(x) | isnotfinite(y)) && return dvi_fpdb_db_nonfinite(x,y)
     return DoubleFloat{T}(dvi_fpdd_dd(x, HILO(y)))
 end
 
 @inline function add_fpdb_db_nonfinite(x::T, y::DoubleFloat{T}) where {T<:IEEEFloat}
     z = x + HI(y)
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function sub_fpdb_db_nonfinite(x::T, y::DoubleFloat{T}) where {T<:IEEEFloat}
     z = x - HI(y)
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function mul_fpdb_db_nonfinite(x::T, y::DoubleFloat{T}) where {T<:IEEEFloat}
     z = x * HI(y)
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
 
 @inline function dvi_fpdb_db_nonfinite(x::T, y::DoubleFloat{T}) where {T<:IEEEFloat}
     z = x / HI(y)
-    return DoubleFloat{T}(z, T(NaN))
+    return DoubleFloat{T}(z, z)
 end
diff --git a/src/type/predicates.jl b/src/type/predicates.jl
@@ -1,4 +1,4 @@
-""""
+"""
 __Predicates__ are functions that ask "yes or no" questions of their argument[s].
 You can ask of a number "Is this zero?" or "Is this one?" and these predicates
 (`iszero`, `isone`) will work as expected with almost all numerical types.
@@ -87,10 +87,10 @@ isinf(x::DoubleFloat{T}) where {T<:IEEEFloat} =
 """
     isposinf(x)
 
-Tests whether a number positive and infinite.
+Tests whether a number is positive and infinite.
 """
 isposinf(x::DoubleFloat{T}) where {T<:IEEEFloat} =
-    isinf(HI(x))
+    isinf(HI(x)) && !signbit(HI(x))
 
 """
     isneginf(x)

diff --git a/src/type/specialvalues.jl b/src/type/specialvalues.jl
@@ -2,9 +2,9 @@ zero(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(zero(T), zero
 one(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(one(T), zero(T))
 
 nan(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(T(NaN), T(NaN))
-inf(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(T(Inf), T(NaN))
-posinf(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(T(Inf), T(NaN))
-neginf(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(T(-Inf), T(NaN))
+inf(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(T(Inf), T(Inf))
+posinf(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(T(Inf), T(Inf))
+neginf(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(T(-Inf), T(-Inf))
 
 typemax(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(typemax(T))
 typemin(::Type{DoubleFloat{T}}) where {T<:IEEEFloat} = DoubleFloat{T}(typemin(T))

diff --git a/test/specialvalues.jl b/test/specialvalues.jl
@@ -13,16 +13,19 @@ end
     @test isinf(T(Inf)) == isinf(inf(T))
     @test isinf(T(Inf)) == isposinf(posinf(T))
     @test isinf(T(-Inf)) == isneginf(neginf(T))
+    @test isinf(T(-Inf)) == !isposinf(neginf(T))
     @test isnan(T(NaN)) == isnan(nan(T))
 
     @test isinf(T(Inf32)) == isinf(inf(T))
     @test isinf(T(Inf32)) == isposinf(posinf(T))
     @test isinf(T(-Inf32)) == isneginf(neginf(T))
+    @test isinf(T(-Inf32)) == !isposinf(neginf(T))
     @test isnan(T(NaN32)) == isnan(nan(T))
 
     @test isinf(T(Inf16)) == isinf(inf(T))
     @test isinf(T(Inf16)) == isposinf(posinf(T))
     @test isinf(T(-Inf16)) == isneginf(neginf(T))
+    @test isinf(T(-Inf16)) == !isposinf(neginf(T))
     @test isnan(T(NaN16)) == isnan(nan(T))
 end
 
@@ -76,6 +79,18 @@ end
     @test isnan(asech(T(NaN)))
     @test isnan(acoth(T(NaN)))
 end
+
+@testset "Inf and NaN function values $T" for T in (Double16, Double32, Double64)
+    @test isnan(sqrt(T(0)/T(0)))
+    @test isinf(sqrt(T(Inf)))
+    @test isinf(exp(T(Inf)))
+    @test isinf(log(T(Inf)))
+    @test isinf(log2(T(Inf)))
+    @test isinf(log10(T(Inf)))
+    @test isinf(log1p(T(Inf)))
+    @test isinf(T(Inf)*T(1))
+    @test isinf(T(1)*T(Inf))
+end
 
 @testset "floatmin2 $T" for T in (Double16, Double32, Double64)
     trueval = (twopar = 2one(T); twopar^trunc(Integer,log(floatmin(T)/eps(T))/log(twopar)/twopar))