Skip to content

Commit

Permalink
1.0-dev: fix interval(::Interval) and some ambiguities (#554)
Browse files Browse the repository at this point in the history
* Fix `interval(::Interval)` and some ambiguities

* Remove superfluous`promote_rule`

* Prevent promotion between `Real` and `Interval`
  • Loading branch information
OlivierHnt authored Feb 14, 2023
1 parent 7892c35 commit 63a14a3
Show file tree
Hide file tree
Showing 4 changed files with 54 additions and 31 deletions.
46 changes: 21 additions & 25 deletions src/intervals/arithmetic/basic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -8,17 +8,9 @@

+(a::Interval) = a # Not in the IEEE standard

"""
-(a::Interval)
Implement the `neg` function of the IEEE Std 1788-2015 (Table 9.1).
"""
-(a::F) where {F<:Interval} = F(-sup(a), -inf(a))


"""
+(a::Interval, b::Real)
+(a::Real, a::Interval)
+(a::Real, b::Interval)
+(a::Interval, b::Interval)
Implement the `add` function of the IEEE Std 1788-2015 (Table 9.1).
Expand All @@ -34,10 +26,18 @@ function +(a::F, b::F) where {F<:Interval}
(isempty(a) || isempty(b)) && return emptyinterval(F)
return @round(F, inf(a) + inf(b), sup(a) + sup(b))
end
+(a::Interval, b::Interval) = +(promote(a, b)...)

"""
-(a::Interval)
Implement the `neg` function of the IEEE Std 1788-2015 (Table 9.1).
"""
-(a::F) where {F<:Interval} = F(-sup(a), -inf(a))

"""
-(a::Interval, b::Real)
-(a::Real, a::Interval)
-(a::Real, b::Interval)
-(a::Interval, b::Interval)
Implement the `sub` function of the IEEE Std 1788-2015 (Table 9.1).
Expand All @@ -46,19 +46,18 @@ function -(a::F, b::T) where {T<:Real, F<:Interval{T}}
isempty(a) && return emptyinterval(F)
return @round(F, inf(a) - b, sup(a) - b)
end

function -(b::T, a::F) where {T, F<:Interval{T}}
isempty(a) && return emptyinterval(F)
return @round(F, b - sup(a), b - inf(a))
end
-(a::F, b::Real) where {F<:Interval} = a - F(b)
-(a::Real, b::F) where {F<:Interval} = F(a) - b

function -(a::F, b::F) where {F<:Interval}
(isempty(a) || isempty(b)) && return emptyinterval(F)
return @round(F, inf(a) - sup(b), sup(a) - inf(b))
end

-(a::F, b::Real) where {F<:Interval} = a - F(b)
-(a::Real, b::F) where {F<:Interval} = F(a) - b
-(a::Interval, b::Interval) = -(promote(a, b)...)

"""
scale(α, a::Interval)
Expand All @@ -71,7 +70,7 @@ For efficiency, does not check that the constant is positive.

"""
*(a::Interval, b::Real)
*(a::Real, a::Interval)
*(a::Real, b::Interval)
*(a::Interval, b::Interval)
Implement the `mul` function of the IEEE Std 1788-2015 (Table 9.1).
Expand All @@ -88,37 +87,33 @@ function *(x::T, a::F) where {T<:Real, F<:Interval{T}}
return @round(F, sup(a)*x, inf(a)*x)
end
end

*(x::T, a::F) where {T<:Real, S, F<:Interval{S}} = Interval{S}(x) * a
*(a::F, x::T) where {T<:Real, S, F<:Interval{S}} = x*a

function *(a::F, b::F) where {F<:Interval}
(isempty(a) || isempty(b)) && return emptyinterval(F)

(isthinzero(a) || isthinzero(b)) && return zero(F)

(isbounded(a) && isbounded(b)) && return mult(*, a, b)

return mult((x, y, r) -> unbounded_mult(F, x, y, r), a, b)
end

*(a::Interval, b::Interval) = *(promote(a, b)...)

# Helper functions for multiplication
function unbounded_mult(::Type{F}, x::T, y::T, r::RoundingMode) where {T, F<:Interval{T}}
iszero(x) && return sign(y)*zero_times_infinity(T)
iszero(y) && return sign(x)*zero_times_infinity(T)
iszero(x) && return sign(y) * zero_times_infinity(T)
iszero(y) && return sign(x) * zero_times_infinity(T)
return *(x, y, r)
end

function mult(op, a::F, b::F) where {T, F<:Interval{T}}
if inf(b) >= zero(T)
inf(a) >= zero(T) && return @round(F, op(inf(a), inf(b)), op(sup(a), sup(b)))
sup(a) <= zero(T) && return @round(F, op(inf(a), sup(b)), op(sup(a), inf(b)))
return @round(F, inf(a)*sup(b), sup(a)*sup(b)) # when zero(T) ∈ a
return @round(F, inf(a)*sup(b), sup(a)*sup(b)) # zero(T) ∈ a
elseif sup(b) <= zero(T)
inf(a) >= zero(T) && return @round(F, op(sup(a), inf(b)), op(inf(a), sup(b)))
sup(a) <= zero(T) && return @round(F, op(sup(a), sup(b)), op(inf(a), inf(b)))
return @round(F, sup(a)*inf(b), inf(a)*inf(b)) # when zero(T) ∈ a
return @round(F, sup(a)*inf(b), inf(a)*inf(b)) # zero(T) ∈ a
else
inf(a) > zero(T) && return @round(F, op(sup(a), inf(b)), op(sup(a), sup(b)))
sup(a) < zero(T) && return @round(F, op(inf(a), sup(b)), op(inf(a), inf(b)))
Expand All @@ -129,7 +124,7 @@ end

"""
/(a::Interval, b::Real)
/(a::Real, a::Interval)
/(a::Real, b::Interval)
/(a::Interval, b::Interval)
Implement the `div` function of the IEEE Std 1788-2015 (Table 9.1).
Expand Down Expand Up @@ -182,6 +177,7 @@ function /(a::F, b::F) where {T, F<:Interval{T}}
end
end
end
/(a::Interval, b::Interval) = /(promote(a, b)...)

"""
inv(a::Interval)
Expand Down
15 changes: 9 additions & 6 deletions src/intervals/construction.jl
Original file line number Diff line number Diff line change
Expand Up @@ -109,6 +109,11 @@ Interval(x::Irrational) = Interval{default_bound()}(x)
return :(return $res) # Set body of the function to return the precomputed result
end

# promotion

Base.promote_rule(::Type{Interval{T}}, ::Type{Interval{S}}) where {T,S} =
Interval{promote_type(T, S)}

"""
interval(a, b)
Expand All @@ -117,15 +122,13 @@ If so, then an `Interval(a, b)` object is returned;
if not, a warning is printed and the empty interval is returned.
"""
function interval(a::T, b::S) where {T<:Real, S<:Real}
if !is_valid_interval(a, b)
@warn "Invalid input, empty interval is returned"
return emptyinterval(promote_type(T, S))
end

return Interval(a, b)
is_valid_interval(a, b) && return Interval(a, b)
@warn "Invalid input, empty interval is returned"
return emptyinterval(promote_type(T, S))
end

interval(a::Real) = interval(a, a)
interval(a::Interval) = interval(inf(a), sup(a)) # Check the validity of the interval

const checked_interval = interval

Expand Down
18 changes: 18 additions & 0 deletions test/interval_tests/construction.jl
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,9 @@ using Test
@test big(ℯ) in Interval{Float32}(0, ℯ)
@test big(π) in Interval{Float32}(π, 4)

@test interval(Interval(pi)) Interval(pi)
@test interval(Interval(NaN, -Inf)) emptyinterval()

# a < Inf and b > -Inf
@test @interval("1e300") Interval(9.999999999999999e299, 1.0e300)
@test @interval("-1e307") Interval(-1.0000000000000001e307, -1.0e307)
Expand Down Expand Up @@ -165,6 +168,16 @@ end
@test convert(Interval{BigFloat}, x) === x
end

@testset "Promotion between intervals" begin
x = Interval{Float64}(π)
y = Interval{BigFloat}(π)
x_, y_ = promote(x, y)

@test promote_type(typeof(x), typeof(y)) == Interval{BigFloat}
@test bounds(x_) == (BigFloat(inf(x), RoundDown), BigFloat(sup(x), RoundUp))
@test y_ y
end

@testset "Typed intervals" begin
@test typeof(@interval Float64 1 2) == Interval{Float64}
@test typeof(@interval 1 2) == Interval{Float64}
Expand All @@ -186,6 +199,11 @@ end

a = convert(Interval{Float64}, @biginterval(3, 4))
@test typeof(a) == Interval{Float64}

pi64, pi32 = Interval{Float64}(pi), Interval{Float32}(pi)
x, y = promote(pi64, pi32)
@test x pi64
@test y Interval{Float64}(pi32)
end

@testset "Interval{T} constructor" begin
Expand Down
6 changes: 6 additions & 0 deletions test/interval_tests/numeric.jl
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,12 @@ end
@test a + b Interval(+(a.lo, b.lo, RoundDown), +(a.hi, b.hi, RoundUp))
@test -a Interval(-a.hi, -a.lo)
@test a - b Interval(-(a.lo, b.hi, RoundDown), -(a.hi, b.lo, RoundUp))
for f in (:+, :-, :*, :/)
@eval begin
@test $f(Interval{Float64}(pi), Interval{Float32}(pi))
$f(Interval{Float64}(pi), Interval{Float64}(Interval{Float32}(pi)))
end
end
@test Interval(1//4,1//2) + Interval(2//3) Interval(11//12, 7//6)
@test_broken Interval(1//4,1//2) - Interval(2//3) Interval(-5//12, -1//6)

Expand Down

0 comments on commit 63a14a3

Please sign in to comment.