Help with IntervalArithmetic

intervals

#1

I am using the IntervalArithmetic packaged and I need to convert a vector and a matrix of floats to a vector and matrix of intervals. Is there a simple way to do this?

More specifically, let’s say I have

v = [1; 2.1]
A = [1 3.2; 2 3]

and I want to transform v and A in a interval vector and an interval matrix. I can do this with the commands

vi = [Interval(v[1]); Interval(v[2])]
Ai = [Interval(A[1,1]) Interval(A[1,2]); Interval(A[2,1]) Interval(A[2,2])]

but this is not convenient for vectors and matrices of arbitrary sizes. Is there a simpler way to do this? It would be nice to be able to do

vi = Interval(v)
Ai = Interval(A)

like in INTLAB.

A second question is, how do I get an interval that is a tight enclosure to a float. If I have, for example x = 2.1 and I type Interval(x) I get the interval

[2.1, 2.10001]

How can I get the interval [x, x] given that x is already a floating point number defined in my code.

Thanks


#2

To answer your first question, you can use the map function.

vi = map(Interval, v)
Ai = map(Interval, A)

I’m not sure about your second question.


#3

Thanks!


#4

More tersely:

julia> Interval.(v)
2-element Array{IntervalArithmetic.Interval{Float64},1}:
 [1, 1]            
     [2.1, 2.10001]

julia> Interval.(A)
2×2 Array{IntervalArithmetic.Interval{Float64},2}:
 [1, 1]      [3.2, 3.20001]
 [2, 2]  [3, 3]

#5

Hi, IntervalArithmetic.jl author here.

If I have, for example x = 2.1 and I type Interval(x) I get the interval
[2.1, 2.10001]

Actually you don’t; this is just the output indicating that the result is contained inside that interval. (I believe that Intlab does the same.) Interval(x) returns the “thin” interval [x, x] (i.e. a single point). You can change the output format to see this:

julia 0.6> @format full
6

julia 0.6> x = 2.1
2.1

julia 0.6> Interval(x)
Interval(2.1, 2.1)

If you want a very tiny interval around x that contains the true real number 2.1, you can do one of the following two things:

julia 0.6> x..x
Interval(2.0999999999999996, 2.1000000000000005)

julia 0.6> convert(Interval, x)
Interval(2.0999999999999996, 2.1000000000000005)

(Note that this output has changed slightly in the latest releases of IntervalArithmetic.jl, compared with the previous releases.)

Just a reminder that when you write x=2.1, the result is not the exact real number 2.1, but is rather
the following exact decimal number:

julia 0.6> big(2.1)
2.100000000000000088817841970012523233890533447265625000000000000000000000000000

This is due to how floating-point arithmetic works.


#6

Thank you for the reply! There is still one thing I want to make sure I understand correctly.

I am aware that when I type x = 2.1, for example, x is not necessarily the number I typed. It would be the number I typed only if the number I typed were exactly representable as a floating point number. Please correct me if I am wrong!

In my case I have floating point numbers that I computed previously and I want to make then into intervals, so Interval(x) is exactly what I need, that is, I need the thin interval containing the floating point number x.

If I want an interval containing 2.1, I understand that I can use

y = 2.1…2.1

and this will make sure 2.1 in within my interval.

Now my question is, if I have a floating point number, x = 2.1 for example, then what exactly

y = x…x

does? How does it guarantee that 2.1 is in this interval? Does it take the interval [x1, x2], where x1 is the largest representable number smaller than x and x2 is the smallest representable number larger than x? If this is the case then I can see that 2.1 would be contained in x…x.

However, if this is the case, the width of the interval I obtain with 2.1…2.1 could be larger than needed when the number I type is exactly representable. Is there a way to type a number a and get the thin interval [a, a] if a is exactly representable and get an interval containing a if a is not exactly representable?

Thanks


#7

I am aware that when I type x = 2.1, for example, x is not necessarily the number I typed. It would be the number I typed only if the number I typed were exactly representable as a floating point number. Please correct me if I am wrong!

That is exactly correct. For example, x = 2.5 is exactly representable.

If I want an interval containing 2.1, I understand that I can use

y = 2.1…2.1

and this will make sure 2.1 in within my interval.

Correct.

Now my question is, if I have a floating point number, x = 2.1 for example, then what exactly

y = x…x

does? How does it guarantee that 2.1 is in this interval? Does it take the interval [x1, x2], where x1 is the largest representable number smaller than x and x2 is the smallest representable number larger than x? If this is the case then I can see that 2.1 would be contained in x…x.

This is almost correct. What it does is form an Interval using the following two numbers:

julia 0.6> (prevfloat(x), nextfloat(x))
(2.0999999999999996, 2.1000000000000005)

That is, one below and one above the floating-point number x.

In previous versions of the package, it tried to be more clever, by e.g. converting x to a rational first:

julia 0.6> rationalize(x)
21//10

julia 0.6> convert(Interval{Float64}, rationalize(x))
Interval(2.0999999999999996, 2.1)

This literally does do what you suggested; however, it is much slower. The version using prevfloat only loses one “ulp”, i.e. it is one floating-point number wider.

However, if this is the case, the width of the interval I obtain with 2.1…2.1 could be larger than needed when the number I type is exactly representable. Is there a way to type a number a and get the thin interval [a, a] if a is exactly representable and get an interval containing a if a is not exactly representable?

You can do the following, but it is quite slow:

julia 0.6> x = 2.1; s = string(x)
"2.1"

julia 0.6> Interval(parse(Float64, s, RoundDown), parse(Float64, s, RoundUp))
Interval(2.0999999999999996, 2.1)

julia 0.6> x = 2.5; s = string(x)
"2.5"

julia 0.6> Interval(parse(Float64, s, RoundDown), parse(Float64, s, RoundUp))
Interval(2.5, 2.5)

#8

Thank you! Now I understand how to set the intervals.

I have two more questions:

Will Interval(pi) and Interval(e) give me intervals that are guaranteed to contain the number pi and the number e? Or do I have to type in something explicit such as 3.141592…3.141593 to get an interval containing pi?

If I multiply an interval by a float or add a float to an interval, for example 2.1*x or x + 2.1, where x is an interval, do I get intervals that are guaranteed to contain the correct results? I know that this is not true for more complicated expressions, and in those cases I can use the macro @intervals( ).


#9

Great questions! If you don’t mind, I may use them as FAQs for future versions of the package.

Will Interval(pi) and Interval(e) give me intervals that are guaranteed to contain the number pi and the number e? Or do I have to type in something explicit such as 3.141592…3.141593 to get an interval containing pi?

Since pi and e are floating-point numbers, just using Interval again gives you thin intervals:

julia 0.6> Interval(pi)
Interval(3.141592653589793, 3.141592653589793)

You can do

julia 0.6> pi..pi
Interval(3.1415926535897927, 3.1415926535897936)

or

julia 0.6> @interval(pi)
Interval(3.141592653589793, 3.1415926535897936)

If I multiply an interval by a float or add a float to an interval, for example 2.1*x or x + 2.1, where x is an interval, do I get intervals that are guaranteed to contain the correct results? I know that this is not true for more complicated expressions, and in those cases I can use the macro @intervals( ).

Yes, you do get correct results: the 2.1 is automatically converted. This is also true for more complicated expressions, provided there is an interval that forces the conversion.

For example,

julia 0.6> (2.1..2.1)*3.1 + (4.1..4.1)*5.1
Interval(27.419999999999984, 27.420000000000012)

julia 0.6> @interval(2.1*3.1 + 4.1*5.1)
Interval(27.419999999999984, 27.420000000000012)

julia 0.6> 2.1*3.1 + (4.1..4.1)*5.1
Interval(27.419999999999984, 27.42000000000001)

The first two are correct. In the last one, the multiplication 2.1*3.1 is done without any guarantees, since that is evaluated using standard floating-point arithmetic before the intervals come into the story.


#10

I should say that the documentation for the package is definitely lacking and not precise enough as regards the answers to your questions. If you have time and energy, help with the docs (via a Pull Request) would be very welcome.


#11

Fell free the use the questions as FAQs.

I am still not sure about pi for example. If I use x = pi I understand that Julia gives me a floating point approximation of pi, and if I type pi…pi I get an interval that contains this numerical approximation of pi. However I will not have any guarantees that the actual number pi is in this interval, since taking [previous float, next float] may not be enough to compensate for the numerical approximation made by Julia. I am correct?

I am working on projects on computer assisted proofs, so I need to have guarantees that the value of pi is in my interval.


#12

I am kind of new to GitHub, but I will take a look if I can manage to pull the code add the questions and answers to the documentation.


#13

Fell free the use the questions as FAQs.

Thanks!

I am still not sure about pi for example. If I use x = pi I understand that Julia gives me a floating point approximation of pi, and if I type pi…pi I get an interval that contains this numerical approximation of pi. However I will not have any guarantees that the actual number pi is in this interval, since taking [previous float, next float] may not be enough to compensate for the numerical approximation made by Julia. I am correct?

Actually it is enough, since when you use pi in a standard arithmetic expression in Julia, the nearest floating point number to the true value of π is used; so when you do prevfloat and nextfloat it is guaranteed to include that true value.

You can also do e.g. const PI = @interval(pi) to get the tightest Float64 interval approximation, or @biginterval(pi) to get a very tight Interval using BigFloats. Unfortunately BigFloats are very slow, though.

I am working on projects on computer assisted proofs, so I need to have guarantees that the value of pi is in my interval.

I would be very interested in hearing more details.


#15

You can just click on one of the files here, then click on the pencil icon, and you can edit the file right there in the GitHub interface. The docs use Markdown syntax.