Fix overflow problems with round and ceil for Normed

This also simplifies the test for rounding functions.

Author: kimikage

src/normed.jl

@@ -235,26 +235,35 @@ abs(x::Normed) = x
 /(x::T, y::T) where {T <: Normed} = convert(T,convert(floattype(T), x)/convert(floattype(T), y))
 
 # Functions
-trunc(x::T) where {T <: Normed} = T(div(reinterpret(x), rawone(T))*rawone(T),0)
-floor(x::T) where {T <: Normed} = trunc(x)
-function round(x::Normed{T,f}) where {T,f}
-    mask = convert(T, 1<<(f-1))
-    y = trunc(x)
-    return convert(T, reinterpret(x)-reinterpret(y)) & mask>0 ?
-            Normed{T,f}(y+oneunit(Normed{T,f})) : y
+trunc(x::N) where {N <: Normed} = floor(x)
+floor(x::N) where {N <: Normed} = reinterpret(N, x.i - x.i % rawone(N))
+function ceil(x::Normed{T,f}) where {T, f}
+    f == 1 && return x
+    if typemax(T) % rawone(x) != 0
+        upper = typemax(T) - typemax(T) % rawone(x)
+        x.i > upper && throw_converterror(Normed{T,f}, ceil(T, typemax(x)))
+    end
+    r = x.i % rawone(x)
+    reinterpret(Normed{T,f}, x.i - r + (r > 0 ? rawone(x) : zero(T)))
 end
-function ceil(x::Normed{T,f}) where {T,f}
-    k = bitwidth(T)-f
-    mask = (typemax(T)<<k)>>k
-    y = trunc(x)
-    return convert(T, reinterpret(x)-reinterpret(y)) & (mask)>0 ?
-            Normed{T,f}(y+oneunit(Normed{T,f})) : y
+function round(x::Normed{T,f}) where {T, f}
+    r = x.i % rawone(x)
+    q = rawone(x) - r
+    reinterpret(Normed{T,f}, r > q ? x.i + q : x.i - r)
 end
 
-trunc(::Type{T}, x::Normed) where {T <: Integer} = convert(T, div(reinterpret(x), rawone(x)))
-round(::Type{T}, x::Normed) where {T <: Integer} = round(T, reinterpret(x)/rawone(x))
-floor(::Type{T}, x::Normed) where {T <: Integer} = trunc(T, x)
- ceil(::Type{T}, x::Normed) where {T <: Integer} =  ceil(T, reinterpret(x)/rawone(x))
+trunc(::Type{Ti}, x::Normed) where {Ti <: Integer} = floor(Ti, x)
+function floor(::Type{Ti}, x::Normed) where {Ti <: Integer}
+    convert(Ti, reinterpret(x) ÷ rawone(x))
+end
+function ceil(::Type{Ti}, x::Normed) where {Ti <: Integer}
+    d, r = divrem(x.i, rawone(x))
+    convert(Ti, r > 0 ? d + oneunit(rawtype(x)) : d)
+end
+function round(::Type{Ti}, x::Normed) where {Ti <: Integer}
+    d, r = divrem(x.i, rawone(x))
+    convert(Ti, r > (rawone(x) >> 0x1) ? d + oneunit(rawtype(x)) : d)
+end
 
 isfinite(x::Normed) = true
 isnan(x::Normed) = false

test/normed.jl

@@ -233,41 +233,29 @@ end
     end
 end
 
-function testtrunc(inc::T) where {T}
-    incf = convert(Float64, inc)
-    tm = reinterpret(typemax(T))/reinterpret(one(T))
-    local x = zero(T)
-    for i = 0 : min(1e6, reinterpret(typemax(T))-1)
-        xf = incf*i
-        try
-            @test typeof(trunc(x)) == T
-            @test trunc(x) == trunc(xf)
-            @test typeof(round(x)) == T
-            @test round(x) == round(xf)
-            cxf = ceil(xf)
-            if cxf < tm
-                @test typeof(ceil(x)) == T
-                @test ceil(x) == ceil(xf)
-            end
-            @test typeof(floor(x)) == T
-            @test floor(x) == floor(xf)
-            @test trunc(Int,x) == trunc(Int,xf)
-            @test round(Int,x) == round(Int,xf)
-            @test floor(Int,x) == floor(Int,xf)
-            if cxf < tm
-                @test ceil(Int,x) == ceil(Int,xf)
-            end
-        catch err
-            println("Failed on x = ", x, ", xf = ", xf)
-            rethrow(err)
+@testset "rounding" begin
+    for T in (UInt8, UInt16, UInt32, UInt64)
+        rs = vcat([ oneunit(T) << b - oneunit(T) << 1 for b = 1:bitwidth(T)],
+                  [ oneunit(T) << b - oneunit(T)      for b = 1:bitwidth(T)],
+                  [ oneunit(T) << b                   for b = 2:bitwidth(T)-1])
+        @testset "rounding Normed{$T,$f}" for f = 1:bitwidth(T)
+            N = Normed{T,f}
+            xs = (reinterpret(N, r) for r in rs)
+            @test all(x -> trunc(x) == trunc(float(x)), xs)
+            @test all(x -> floor(x) == floor(float(x)), xs)
+            # force `Normed` comparison avoiding rounding errors
+            @test all(x -> ceil(float(x)) > typemax(N) || ceil(x) == N(ceil(float(x))), xs)
+            @test all(x -> round(x) == round(float(x)), xs)
+            @test all(x -> trunc(UInt64, x) === trunc(UInt64, float(x)), xs)
+            @test all(x -> floor(UInt64, x) === floor(UInt64, float(x)), xs)
+            @test all(x ->  ceil(UInt64, x) ===  ceil(UInt64, float(x)), xs)
+            @test all(x -> round(UInt64, x) === round(UInt64, float(x)), xs)
         end
-        x = convert(T, x+inc)
     end
-end
-
-@testset "trunc" begin
-    for T in (FixedPointNumbers.UF..., UF2...)
-        testtrunc(eps(T))
+    @testset "rounding Normed{Int16,$f} with overflow" for f in filter(x->!ispow2(x), 1:16)
+        N = Normed{UInt16,f}
+        @test_throws ArgumentError ceil(typemax(N))
+        @test_throws ArgumentError ceil(floor(typemax(N)) + eps(N))
     end
 end