From 104a83622d9d2f2aa34230a5f45a199836aec30e Mon Sep 17 00:00:00 2001 From: David Sanders Date: Sat, 23 Feb 2019 08:40:41 -0600 Subject: [PATCH 1/7] Fix some type stability issues --- src/intervals/arithmetic.jl | 20 ++++++++++---------- 1 file changed, 10 insertions(+), 10 deletions(-) diff --git a/src/intervals/arithmetic.jl b/src/intervals/arithmetic.jl index 00dc2d092..f0fc61c08 100644 --- a/src/intervals/arithmetic.jl +++ b/src/intervals/arithmetic.jl @@ -162,8 +162,8 @@ function inv(a::Interval{T}) where T<:Real isempty(a) && return emptyinterval(a) if zero(T) ∈ a - a.lo < zero(T) == a.hi && return @round(-Inf, inv(a.lo)) - a.lo == zero(T) < a.hi && return @round(inv(a.hi), Inf) + a.lo < zero(T) == a.hi && return @round(T(-Inf), inv(a.lo)) + a.lo == zero(T) < a.hi && return @round(inv(a.hi), T(Inf)) a.lo < zero(T) < a.hi && return entireinterval(T) a == zero(a) && return emptyinterval(T) end @@ -195,14 +195,14 @@ function /(a::Interval{T}, b::Interval{T}) where T<:Real if iszero(b.lo) - a.lo >= zero(T) && return @round(a.lo/b.hi, Inf) - a.hi <= zero(T) && return @round(-Inf, a.hi/b.hi) + a.lo >= zero(T) && return @round(a.lo/b.hi, T(Inf)) + a.hi <= zero(T) && return @round(T(-Inf), a.hi/b.hi) return entireinterval(S) elseif iszero(b.hi) - a.lo >= zero(T) && return @round(-Inf, a.lo/b.lo) - a.hi <= zero(T) && return @round(a.hi/b.lo, Inf) + a.lo >= zero(T) && return @round(T(-Inf), a.lo/b.lo) + a.hi <= zero(T) && return @round(a.hi/b.lo, T(Inf)) return entireinterval(S) else @@ -218,10 +218,10 @@ function extended_div(a::Interval{T}, b::Interval{T}) where T<:Real S = typeof(a.lo / b.lo) if 0 < b.hi && 0 > b.lo && 0 ∉ a if a.hi < 0 - return (Interval(-Inf, a.hi / b.hi), Interval(a.hi / b.lo, Inf)) + return (Interval(T(-Inf), a.hi / b.hi), Interval(a.hi / b.lo, T(Inf))) elseif a.lo > 0 - return (Interval(-Inf, a.lo / b.lo), Interval(a.lo / b.hi, Inf)) + return (Interval(T(-Inf), a.lo / b.lo), Interval(a.lo / b.hi, T(Inf))) end elseif 0 ∈ a && 0 ∈ b @@ -432,13 +432,13 @@ function mid(a::Interval{T}) where T a.lo == -∞ && return nextfloat(a.lo) a.hi == +∞ && return prevfloat(a.hi) - midpoint = 0.5 * (a.lo + a.hi) + midpoint = (a.lo + a.hi) / 2 isfinite(midpoint) && return midpoint #= Fallback in case of overflow: a.hi + a.lo == +∞ or a.hi + a.lo == -∞. This case can not be the default one as it does not pass several IEEE1788-2015 tests for small floats. =# - return 0.5 * a.lo + 0.5 * a.hi + return a.lo / 2 + a.hi / 2 end mid(a::Interval{Rational{T}}) where T = (1//2) * (a.lo + a.hi) From ef08499b134a46c245ad2e517067f7dd9753622c Mon Sep 17 00:00:00 2001 From: David Sanders Date: Sat, 23 Feb 2019 11:30:11 -0600 Subject: [PATCH 2/7] Fix escaping in @round --- src/intervals/rounding_macros.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/intervals/rounding_macros.jl b/src/intervals/rounding_macros.jl index 751a934de..ddd079740 100644 --- a/src/intervals/rounding_macros.jl +++ b/src/intervals/rounding_macros.jl @@ -22,7 +22,7 @@ function round_expr(ex::Expr, rounding_mode::RoundingMode) return :( $(esc(op))( $(esc(ex.args[2])), $(esc(ex.args[3])), $rounding_mode) ) else # unary operator - return :( $op($(esc(ex.args[2])), $rounding_mode ) ) + return :( $(esc(op))($(esc(ex.args[2])), $rounding_mode ) ) end else return ex From d97d79b367818513e4c9e91c8b4fd6e5e97469f3 Mon Sep 17 00:00:00 2001 From: David Sanders Date: Sun, 24 Feb 2019 01:40:33 -0600 Subject: [PATCH 3/7] Make mid type stable --- src/intervals/arithmetic.jl | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/src/intervals/arithmetic.jl b/src/intervals/arithmetic.jl index f0fc61c08..36ac59485 100644 --- a/src/intervals/arithmetic.jl +++ b/src/intervals/arithmetic.jl @@ -394,7 +394,7 @@ The default is the true midpoint at `α = 0.5`. Assumes 0 ≤ α ≤ 1. -Warning: if the parameter `α = 0.5` is explicitely set, the behavior differs +Warning: if the parameter `α = 0.5` is explicitly set, the behavior differs from the default case if the provided `Interval` is not finite, since when `α` is provided `mid` simply replaces `+∞` (respectively `-∞`) by `prevfloat(+∞)` (respecively `nextfloat(-∞)`) for the computation of the intermediate point. @@ -406,13 +406,15 @@ function mid(a::Interval{T}, α) where T lo = (a.lo == -∞ ? nextfloat(-∞) : a.lo) hi = (a.hi == +∞ ? prevfloat(+∞) : a.hi) - midpoint = α * (hi - lo) + lo + β = convert(T, α) + + midpoint = β * (hi - lo) + lo isfinite(midpoint) && return midpoint #= Fallback in case of overflow: hi - lo == +∞. This case can not be the default one as it does not pass several IEEE1788-2015 tests for small floats. =# - return (1-α) * lo + α * hi + return (1 - β) * lo + β * hi end """ From 5436e9a037b9656fdd1f1793fa834ad49198d93c Mon Sep 17 00:00:00 2001 From: David Sanders Date: Sun, 24 Feb 2019 01:55:37 -0600 Subject: [PATCH 4/7] Another mid fix --- src/intervals/arithmetic.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/intervals/arithmetic.jl b/src/intervals/arithmetic.jl index 36ac59485..bdacb7f0b 100644 --- a/src/intervals/arithmetic.jl +++ b/src/intervals/arithmetic.jl @@ -403,8 +403,8 @@ function mid(a::Interval{T}, α) where T isempty(a) && return convert(T, NaN) - lo = (a.lo == -∞ ? nextfloat(-∞) : a.lo) - hi = (a.hi == +∞ ? prevfloat(+∞) : a.hi) + lo = (a.lo == -∞ ? nextfloat(T(-∞)) : a.lo) + hi = (a.hi == +∞ ? prevfloat(T(+∞)) : a.hi) β = convert(T, α) From 31cefa4eb6f9bdcb69bfbc85f64e39cd8a2832df Mon Sep 17 00:00:00 2001 From: David Sanders Date: Mon, 4 Mar 2019 19:27:48 -0600 Subject: [PATCH 5/7] Add tests of type stability using @inferred --- test/interval_tests/consistency.jl | 16 ++++++++++++++++ 1 file changed, 16 insertions(+) diff --git a/test/interval_tests/consistency.jl b/test/interval_tests/consistency.jl index 1052e7586..3152696b3 100644 --- a/test/interval_tests/consistency.jl +++ b/test/interval_tests/consistency.jl @@ -388,4 +388,20 @@ setprecision(Interval, Float64) end + @testset "Type stability" begin + for T in (Float32, Float64, BigFloat) + + x = Interval{T}(3, 4) + y = Interval{T}(5, 6) + + for op in (+, -, *, /, atan) + @inferred op(a, b) + end + + for op in (sin, cos, exp, log, tan, abs) + @inferred op(a) + end + end + end + end From 4a2183f33863452227d2d43e9c75c90ac6667193 Mon Sep 17 00:00:00 2001 From: David Sanders Date: Mon, 4 Mar 2019 19:31:31 -0600 Subject: [PATCH 6/7] Add more functions --- test/interval_tests/consistency.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/interval_tests/consistency.jl b/test/interval_tests/consistency.jl index 3152696b3..8cd161d4b 100644 --- a/test/interval_tests/consistency.jl +++ b/test/interval_tests/consistency.jl @@ -398,7 +398,7 @@ setprecision(Interval, Float64) @inferred op(a, b) end - for op in (sin, cos, exp, log, tan, abs) + for op in (sin, cos, exp, log, tan, abs, mid, diam) @inferred op(a) end end From eac765264e1e2fa33db83e3dbdab76d50b210141 Mon Sep 17 00:00:00 2001 From: David Sanders Date: Tue, 5 Mar 2019 13:21:46 -0600 Subject: [PATCH 7/7] Improve tests --- test/interval_tests/consistency.jl | 20 +++++++++++++------- 1 file changed, 13 insertions(+), 7 deletions(-) diff --git a/test/interval_tests/consistency.jl b/test/interval_tests/consistency.jl index 8cd161d4b..15d3c7280 100644 --- a/test/interval_tests/consistency.jl +++ b/test/interval_tests/consistency.jl @@ -391,15 +391,21 @@ setprecision(Interval, Float64) @testset "Type stability" begin for T in (Float32, Float64, BigFloat) - x = Interval{T}(3, 4) - y = Interval{T}(5, 6) + xs = [3..4, 0..4, 0..0, -4..0, -4..4, -Inf..4, 4..Inf, -Inf..Inf] - for op in (+, -, *, /, atan) - @inferred op(a, b) - end + for x in xs + for y in xs + xx = Interval{T}(x) + yy = Interval{T}(y) + + for op in (+, -, *, /, atan) + @inferred op(x, y) + end + end - for op in (sin, cos, exp, log, tan, abs, mid, diam) - @inferred op(a) + for op in (sin, cos, exp, log, tan, abs, mid, diam) + @inferred op(x) + end end end end