From b37bc9af12d0d981428c00b9ac722e82c6b3b187 Mon Sep 17 00:00:00 2001 From: Oscar Smith Date: Sun, 20 Aug 2023 19:15:53 -0400 Subject: [PATCH 1/2] fix UB in `expm1(::Union{Float16, Float32}) `unsafe_trunc(UInt, -1.0)` is UB but worked fine on apple and AMD so we didn't notice??? This has been very broken since 1.7 --- base/special/exp.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/base/special/exp.jl b/base/special/exp.jl index 9cca6f568305f..8e940a4d85ad9 100644 --- a/base/special/exp.jl +++ b/base/special/exp.jl @@ -460,7 +460,7 @@ function expm1(x::Float32) end x = Float64(x) N_float = round(x*Ln2INV(Float64)) - N = unsafe_trunc(UInt64, N_float) + N = unsafe_trunc(Int64, N_float) r = muladd(N_float, Ln2(Float64), x) hi = evalpoly(r, (1.0, .5, 0.16666667546642386, 0.041666183019487026, 0.008332997481506921, 0.0013966479175977883, 0.0002004037059220124)) @@ -477,7 +477,7 @@ function expm1(x::Float16) return Float16(x*evalpoly(x, (1f0, .5f0, 0.16666628f0, 0.04166785f0, 0.008351848f0, 0.0013675707f0))) end N_float = round(x*Ln2INV(Float32)) - N = unsafe_trunc(UInt32, N_float) + N = unsafe_trunc(Int32, N_float) r = muladd(N_float, Ln2(Float32), x) hi = evalpoly(r, (1f0, .5f0, 0.16666667f0, 0.041665863f0, 0.008333111f0, 0.0013981499f0, 0.00019983904f0)) small_part = r*hi From b6aaa5741610992a30263a6e338d63fd79ec226f Mon Sep 17 00:00:00 2001 From: Oscar Smith Date: Sun, 20 Aug 2023 19:52:18 -0400 Subject: [PATCH 2/2] add test --- test/math.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/test/math.jl b/test/math.jl index c6600f79b424c..df5c0d23d37a4 100644 --- a/test/math.jl +++ b/test/math.jl @@ -188,6 +188,7 @@ end @test exp10(x) ≈ exp10(big(x)) @test exp2(x) ≈ exp2(big(x)) @test expm1(x) ≈ expm1(big(x)) + @test expm1(T(-1.1)) ≈ expm1(big(T(-1.1))) @test hypot(x,y) ≈ hypot(big(x),big(y)) @test hypot(x,x,y) ≈ hypot(hypot(big(x),big(x)),big(y)) @test hypot(x,x,y,y) ≈ hypot(hypot(big(x),big(x)),hypot(big(y),big(y)))