A Standards compliant Interval Arithmetic Library

Solving some bugs with the help of ITF1788 Test Suite to make the package a little more close to perfect is what I have worked on this summer. From "I can't solve this error, I have no idea about it" to "Hurrah!, I have finally solved it" it has really been an enriching experience for me. The way I have approached my endeavor is by taking every new error as a challenge and no matter how many hours(or days) it may take to learn about that error, it has to be solved.

using Pkg;

Pkg.add("IntervalArithmetic")

Brief Background

The main focus of interval arithmetic is a simple way to calculate upper and lower endpoints for the range of values of a function in one or more variables. These endpoints are not necessarily the supremum or infimum, since the precise calculation of those values can be difficult or impossible. Treatment is typically limited to real intervals, so quantities of form

ITF1788 stands for Interval Test framework for IEEE 1788, these tests are also written in Julia, and my job is to check what test is not passing in IntervalArithmetic.jl and IntervalContractors.jl and edit the packages to satisfy them. Debugging mathematical code was fun but challenging too. Since a part of my work, in these three months was on decorations and reverse functions, it is good you know about them before reading the main part of the blog.

Decorations

The reverse function of a function f calculates the (interval hull of) its inverse function, f⁻¹. For example, sin_rev(Y::Interval, X::Interval) calculates the (interval hull of) those x ∈ X such that sin(x) ∈ Y. This can also be thought of as an inverse function, calculating X_new := sin⁻¹(Y) ∩ X. The return value is (Y, X_new). Functions such as mul_rev(C::Interval, A::Interval, B::Interval) take three arguments, and correspond to C = A * B; they return (C, A_new, B_new), with A_new and B_new similarly defined to be the corresponding inverse images of the multiplication operator in each component.

I have named my Pull requests according to the name of the files in the test suite.

A. Libieeep1788testsrev.jl

I decided to start with the biggest file in the suite which is libieeep1788testsrev.jl and the first thing I needed to do was add decorations for all the reverse functions in IntervalContractors.jl

Pkg.add(PackageSpec(path="https://github.com/krish8484/IntervalContractors.jl.git", rev="decorated_reverse"))

using IntervalContractors

using IntervalArithmetic

Add decorations for reverse functions

I was a bit confused about how to add decorations for all the functions in IntervalContractors.jl but nothing seemed to work, so that encouraged me to fall back to the standard and then I found this

This solved most of my problems and simply a trivial decoration was provided for these functions. Example,

tan_rev(DecoratedInterval(Interval(-0x1.d02967c31p+53, 0x1.d02967c31p+53), def), DecoratedInterval(Interval(-Inf, 1.5707965), def))

cosh_rev(DecoratedInterval(Interval(1.0, Inf), def), DecoratedInterval(Interval(0.0, Inf), dac))

But the decorations were not the only thing to be implemented, there were some minor bugs in the functions as well(now solved) which are mentioned below

Power_rev when the power is 0 (a = x^0)

There is an issue about this in the repository which was opened last year issue #25 and of course, it was pointed out by the tests as well.

After observing the outputs from other standard-compliant libraries, and thinking about the possible solution, it was finally fixed.

Now, the function returns expected outputs

power_rev(Interval(1.0, 1.0), Interval(1.0, 1.0), 0)[2]

power_rev(Interval(-1.0, 0.0), Interval(5.0, 10.0), 0)[2]

Function mulrevIEEE1788 was not returning the correct outputs

This issue was fixed by using extended_div(a function in IntervalArithmetic) in place of the '/' operator (obviously, not very easy to spot and this required a lot of time)

mul_rev_IEEE1788(Interval(-Inf, -0.1), Interval(-2.1, -0.4))

mul_rev_IEEE1788(Interval(-2.0, 1.1), Interval(0.04, Inf), Interval(0.04, Inf))

Cosh_rev was not returning negative values

Apparently, this function must return negative values as well which it was not and now fixed

cosh_rev(DecoratedInterval(Interval(0x1.8b07551d9f55p+0, 0x1.89bca168970c6p+432), com), DecoratedInterval(Interval(-Inf, 0.0), dac))

cosh_rev(DecoratedInterval(Interval(0x1.8b07551d9f55p+0, 0x1.89bca168970c6p+432), com))

Apart from the above points, a default value was added to each of these functions which were required by the tests This completed libieeep1788testsrev.jl .

All of the above points have been incorporated into one PR(#25).

B. Ieee1788-constructors.jl

Next file to be solved was ieee1788-constructors.jl

Making nai() a decorated interval

According to the standard and the tests(obviously), nai() should be a decorated interval, there was also an open issue(#233) regarding this. Now,

import Pkg; Pkg.add("IntervalArithmetic")

using IntervalArithmetic

@format full

nai()

String Parsing Minor bugs in parsing strings

1)empty in @interval

Now,

@interval("[ empty ]")

2)com was returned as a decoration instead of def, even after it was specified as an argument.

Now,

@decorated("3.56?1_def")

A minor change to deprecate interval_part to interval was also done in the same PR which was an open Issue#212.

After these changes all test from ieee1788-constructors.jl pass and the work is in PR#386. 

Closes Issue#222.

C. Mpfi.jl and Fi.lib.jl

Some new functions were implemented to complete these 2 files from the test suite, the functions are listed below, the only tough part is dealing with these files was to take care of the roundings (Phew!)

1)cot

cot(Interval(0x1.921fb745205b6p+4, 0x1.9ca9d998c8f0cp+4))

2)coth

coth(Interval(-0x3.f03d06503caa2p+92, -0x5.9d7eeea9b9ee0p-264))

There was an open issue(Issue#301) about the roundings of cot and coth from the previous year's JSOC/GSOC which has been solved.

3)cbrt

cbrt(Interval(-0x1856e4be527197p-354, 0xd8p0))

Closed Issue#310

4)csc

csc(Interval(7.0, 7.0))

5)csch

csch(Interval(0x10000000000001p-58, 0x10000000000001p-53))

6)hypot

hypot(Interval(0x1854bfb363dc39p-50, 0x19f625847a5899p-48), Interval(0x1854bfb363dc39p-50, 0x19f625847a5899p-48))

Closed Issue#368

7)sec

sec(Interval(-4.0, -4.0))

8)sech

sech(Interval(0.0, +8.0))

Above work can be found at PR#390

D. C-xsc.jl

The only thing needed to be done in the file is to add a new function which is the nthroot function, it was added and now all tests from this file satisfy...

julia> nthroot(Interval(1024.0, 1024.0), 10)

[2, 2]

julia> nthroot(Interval(0, 81), -4)

[0.333333, ∞]

julia> nthroot(Interval(-81, 81), 4)

[0, 3]

E. Display Interval{Float32} with f0

Pkg.add("IntervalArithmetic")

using IntervalArithmetic

This is a simple display feature to display Float32 with f0 while using Intervals in IntervalArithemtic.jl to solve Issue #387

Interval{Float32}(2, 7)

Above work has been implemented in PR #391

F. Add arithmetic operations of numbers with different types

Earlier, while doing arithmetic operations on Intervals, functions from Base were getting called instead of functions from IntervalArithmetic.jl.

Now, the correct functions are getting called,

julia> @which (3f0..4f0) / 2.1

/(a::Interval{T}, x::S) where {T, S<:AbstractFloat} in IntervalArithmetic at /home/krish8484/.julia/dev/IntervalArithmetic/src/intervals/arithmetic.jl:152

julia> Interval{BigFloat}(2,4)/2

[1, 2]₂₅₆

julia> x = (3f0..4f0) / 2.1

[1.42857f0, 1.90477f0]

julia> @which (3f0..4f0) + 2.1

+(a::Interval{T}, b::S) where {T, S<:AbstractFloat} in IntervalArithmetic at /home/krish8484/.julia/dev/IntervalArithmetic/src/intervals/arithmetic.jl:69

This solves issue #389

Above work has been implemented in PR #392

G. libieeep1788_class.jl

@decorated now responds to setprecision()

julia> typeof(@decorated(1, 2))

DecoratedInterval{Float64}

julia> setprecision(Interval, 256)

256

julia> typeof(@decorated(1, 2))

DecoratedInterval{BigFloat}

@decorated now throws an ArgumentError just like @interval for improper inputs

julia> @decorated(2, -1)

ERROR: ArgumentError: `[2.0, -1.0]` is not a valid interval. Need `a ≤ b` to construct `interval(a, b)`.

Stacktrace:

 [1] interval(::Float64, ::Float64) at C:\Users\ASUS\.julia\dev\IntervalArithmetic\src\intervals\intervals.jl:107

 [2] top-level scope at REPL[2]:1

julia> @decorated("[  -1.0  , 1.0]_fooo") #this was one of ITF1788 Test

ERROR: LoadError: ArgumentError: Cannot process fooo as decoration

Stacktrace:

 [1] parse(::Type{DecoratedInterval{Float64}}, ::String) at C:\Users\ASUS\.julia\dev\IntervalArithmetic\src\parsing.jl:48

 [2] @decorated(::LineNumberNode, ::Module, ::Vararg{Any,N} where N) at C:\Users\ASUS\.julia\dev\IntervalArithmetic\src\decorations\intervals.jl:138

in expression starting at REPL[3]:1

Better string handling for @interval and @decorated

julia> @decorated("[ Nai  ]") == nai()

true

julia> @decorated("[-1.0, +infinity]_def") == DecoratedInterval(Interval(-1.0, Inf), def)

true

julia> @test_throws ArgumentError @interval("[ -inf ,  INF ]_com")

Test Passed

      Thrown: ArgumentError

Error Handling in intervals

Appropriate error is thrown when the exception of that kind occurs

julia> @test_throws IntervalArithmetic.IntvlPartOfNaI interval(@decorated("[nai]"))

Test Passed

      Thrown: IntervalArithmetic.IntvlPartOfNaI

Above work is implemented in PR#398

H. libieeep1788_reduction.jl

Some polishing was to Yash Gupta's Pull Request to add some new functions in the library

julia> vector_sumAbs([1.0, -2.0, 3.0], RoundNearest)

6.0

julia> vector_sum([1.0, 2.0, 3.0], RoundNearest)

6.0

julia> vector_sumSquare([1.0, 2.0, 3.0], RoundNearest)

14.0

julia> vector_dot([1.0, 2.0, 3.0], [1.0, 2.0, 3.0], RoundNearest)

14.0

julia> isnan(vector_dot([1.0, 2.0, -Inf, 4.0], [1.0, 2.0, 0.0, 3.0], RoundNearest))

true

Above work has been implemented in PR #399

I. libieeep1788_bool.jl

Most of the tests from this file were already satisfied. The only thing added was for these functions to deal with nai

Now, all tests from this file pass

julia> strictprecedes(nai(), DecoratedInterval(Interval(3.0, 4.0), trv))

false

julia> isdisjoint(nai(), DecoratedInterval(Interval(3.0, 4.0), def))

false

Above work has been implemented in PR #400

J. libieeep1788_tests_overlap.jl

Most of the work in this file was done by Yash Gupta in PR #266, I have only increased the code coverage and removed the superfluous new file.

This function(overap) gives the relation between 2 intervals

julia> overlap(∅, Interval(1.0, 2.0))

"firstEmpty"

julia> overlap(Interval(1.0, 2.0), Interval(2.0, Inf)) == "meets"

true

julia> overlap(Interval(2.0, 2.0), Interval(0.0, Inf))

"containedBy"

julia> overlap(Interval(1.0, 2.0), Interval(0.0, 2.0))

"finishes"

julia> overlap(Interval(3.0, 4.0), Interval(1.0, 2.0))

"after"

julia> overlap(Interval(1.5, 2.5), Interval(1.0, 2.0))

"overlappedBy"

julia> overlap(Interval(0.0, 3.0), Interval(2.0, 2.0))

"contains"

Concluding Remarks

I have learned a lot this summer, working under David Sanders and trying to debug IntervalArithmetic.jl, I have implemented a lot more than I thought I could and I wish to stay with Julia in the coming years. I hope I get more opportunities in the future to collaborate with mentors for Julia Projects, and in the process be a more prominent part of this wonderful community.