Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Simplified version of mod(x::Interval, y::Real) #525

Merged
merged 12 commits into from
Jul 3, 2023
1 change: 1 addition & 0 deletions src/IntervalArithmetic.jl
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,7 @@ import Base:
in, zero, one, eps, typemin, typemax, abs, abs2, real, min, max,
sqrt, exp, log, sin, cos, tan, cot, inv, cbrt, csc, hypot, sec,
exp2, exp10, log2, log10,
mod,
asin, acos, atan,
sinh, cosh, tanh, coth, csch, sech, asinh, acosh, atanh, sinpi, cospi,
union, intersect, isempty,
Expand Down
18 changes: 18 additions & 0 deletions src/intervals/functions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -373,3 +373,21 @@ function nthroot(a::Interval{T}, n::Integer) where T
b = nthroot(bigequiv(a), n)
return convert(Interval{T}, b)
end

"""
Calculate `x::Interval mod y::Real` where y != zero(y) and y is not inteval`.
"""
function mod(x::Interval, y::Real)
@assert y != zero(y) """mod(x::Interval, y::Real)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

an alternative would be to return an empty interval, I think this would be more conformant with the IEEE standard behavior, e.g.

julia> sqrt(-1.. -1)
∅

julia> log(-2.. -1)
∅

Copy link
Contributor Author

@petvana petvana May 27, 2022

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Or we can define mod(a..b, 0) == 0..0 ⊇ ∅ as Wolfram says
https://www.wolframalpha.com/input?i=limit+x-%3E0%2B+mod%281%2C+x%29

Consequently, for 0 ∈ y, we can define smth like mod(x,y) = union(mod(x,(y.lo..y.lo/3)), mod(x, (y.hi/3)..(y.hi))) because the maximum value needs to by in (y.hi/2)..(y.hi)

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I propose to solve this in some future PR, and move the discussion to #129, to keep track once this PR merged.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Or we can define mod(a..b, 0) == 0..0 ⊇ ∅ as Wolfram says https://www.wolframalpha.com/input?i=limit+x-%3E0%2B+mod%281%2C+x%29

Julia (and Wolfram) return NaN (undefined) for mod(3.2, 0.0), which makes sense due to the division in the algorithm.

Consequently, for 0 ∈ y, we can define smth like mod(x,y) = union(mod(x,(y.lo..y.lo/3)), mod(x, (y.hi/3)..(y.hi))) because the maximum value needs to by in (y.hi/2)..(y.hi)

This idea is interesting, and reminds me a little bit what we have for extended_div. That certainly should be addressed in another PR. Yet, I think for (mathematical) reals, mod(x,y) is a number in the interval [0,y] for y>0, or [y, 0] if y<0, isn't it?

From the range of mod(x,y) for real y, then for a strictly positive interval y, we should return [0, sup(y)], or if it is strictly negative [inf(y), 0]. If 0 ∈ y, we should return the hull of these intervals.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think my last comment corresponds to these results... Or not?

is currently implemented only for a strictly positive or negative divisor y."""
division = x / y
fl = floor(division)
if !isthin(fl)
return y > zero(y) ? Interval(zero(y), y) : Interval(y, zero(y))
else
return y * (division - fl)
end
end

mod(x::Interval, y::Interval) = throw(ArgumentError("mod not defined for interval as divisor `y`"))
mod(x::Real, y::Interval) = throw(ArgumentError("mod not defined for interval as divisor `y`"))
41 changes: 41 additions & 0 deletions test/interval_tests/numeric.jl
Original file line number Diff line number Diff line change
Expand Up @@ -434,3 +434,44 @@ end
@test nthroot(Interval{BigFloat}(-81, -16), -4) == ∅
@test nthroot(Interval{BigFloat}(-81, -16), 1) == Interval{BigFloat}(-81, -16)
end

# approximation used for testing (not to rely on ≈ for intervals)
# ⪆(x, y) = (x ≈ y) && (x ⊇ y)
⪆(x::Interval, y::Interval) = x.lo ≈ y.lo && x.hi ≈ y.hi && x ⊇ y
petvana marked this conversation as resolved.
Show resolved Hide resolved

@testset "`mod`" begin
petvana marked this conversation as resolved.
Show resolved Hide resolved
r = 0.0625
x = r..(1+r)
@test mod(x, 1) == mod(x, 1.0) == 0..1
@test mod(x, 2) == mod(x, 2.0) ⪆ x
@test mod(x, 2.5) ⪆ x
@test mod(x, 0.5) == 0..0.5
@test mod(x, -1) == mod(x, -1.0) == -1..0
@test mod(x, -2) == mod(x, -2.0) ⪆ -2+x
@test mod(x, -2.5) ⪆ -2.5+x
@test mod(x, -0.5) == -0.5..0

x = (-1+r) .. -r
@test mod(x, 1) == mod(x, 1.0) ⪆ 1+x
@test mod(x, 2) == mod(x, 2.0) ⪆ 2+x
@test mod(x, 2.5) ⪆ 2.5+x
@test mod(x, 0.5) == 0..0.5
@test mod(x, -1) == mod(x, -1.0) ⪆ x
@test mod(x, -2) == mod(x, -2.0) ⪆ x
@test mod(x, -2.5) ⪆ x
@test mod(x, -0.5) == -0.5..0

x = -r .. 1-r
@test mod(x, 1) == mod(x, 1.0) == 0..1
@test mod(x, 2) == mod(x, 2.0) == 0..2
@test mod(x, 2.5) == 0..2.5
@test mod(x, 0.5) == 0..0.5
@test mod(x, -1) == mod(x, -1.0) == -1..0
@test mod(x, -2) == mod(x, -2.0) == -2..0
@test mod(x, -2.5) == -2.5..0
@test mod(x, -0.5) == -0.5..0

# TODO - implement mod for two intervals
@test_throws ArgumentError mod(1..2, 1.4..1.5)
@test_throws ArgumentError mod(1.0, 1.4..1.5)
end