High precision interval rounding error

General Usage
,
1

Before everything I need a disclaimer: I do understand that the precision which I set is comically large. However, there is a legitimate reason for what I am doing, but I’m not going to go into it.

The following snippet of code

using IntervalArithmetic
setprecision(8384);
y = @biginterval sqrt(big(2));
x = @biginterval 1;
exponent(diam(y))
exponent(diam(x))

spits out -8383 and then throws ERROR: DomainError with 0.0: "x" must be non-zero and finite.

I assume this is a bug, but I’m not sure where the bug is. It could be in InvervalArithmetic, but it could also be as deep as in MPFR. Does anyone have any insight to this?

2

The problem is that the logical epsilon for 1 is 0 (Int’s are exact). Did you want x = @biginterval big(1.0)?

1 Like
3

It throws the same error.

4

This has nothing to do with intervals as such.
@biginterval 1 contains an interval of 0 width, and then:


julia> exponent(0.0)
ERROR: DomainError with 0.0:
Cannot be subnormal converted to 0.
Stacktrace:
 [1] (::Base.Math.var"#throw2#4")(::Float64) at ./math.jl:784
 [2] exponent(::Float64) at ./math.jl:789
 [3] top-level scope at REPL[7]:100:

It does not make sense to ask for the floating-point exponent of 0.0.

3 Likes
5

In principle intervals with that precision should work fine. There may be some issues with displaying the results; if so, please raise them at https://github.com/JuliaIntervals/IntervalArithmetic.jl/ (and preferably help to solve them :wink: )

2 Likes
6

Oh, I see. Since the number 1 can be represented exactly as a floating point number, the upper and lower bounds of the interval are equal. So then the diameter is exactly 0 and its exponent is not defined.

I assumed that since the floating point number 0 is stored in the memory with some exponent, exponent would return a number. But, I see why this could be problematic in some cases. Thanks!

1 Like
7

0.0 (and -0.0) are both kind of special cased in the floating point standards. If you didn’t, they would be interpreted as 1.0*2^-1024 (except maybe not cause subnormal numbers)

2 Likes