## Interval Arithmetic

### December 21, 2010

Over at his blog, where he has been solving the exercises in SICP, Bill Cruise reached the section of the book that discusses interval arithmetic; if you want an introduction to the topic, you might want to look at http://www.cs.utep.edu/interval-comp/hayes.pdf.

Given an interval [x,y] with x the lower bound and y the upper bound, the basic operations of interval arithmetic are defined as follows:

• [a,b] + [c,d] = [a+c,b+d]
• [a,b] − [c,d] = [ad,bc]
• [a,b] × [c,d] = [min(ac,ad,bc,bd), max(ac,ad,bc,bd)]
• [a,b] ÷ [c,d] = [min(a/c,a/d,b/c,b/d), max(a/c,a/d,b/c,b/d)], where division by an interval containing zero is undefined

Your task is to write a basic library for interval arithmetic; you may do as many of the SICP exercises as you wish. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

### 10 Responses to “Interval Arithmetic”

1. […] Praxis – Interval Arithmetic By Remco Niemeijer In today’s Programming Praxis exercise, our goal is to implement some functions to do interval arithmetic. […]

2. My Haskell solution (see http://bonsaicode.wordpress.com/2010/12/21/programming-praxis-interval-arithmetic/ for a version with comments):

```plus :: (Num a, Num b) => (a, b) -> (a, b) -> (a, b)
plus (a,b) (c,d) = (a+c, b+d)

minus :: (Num a, Num b) => (a, b) -> (b, a) -> (a, b)
minus (a,b) (c,d) = (a-d, b-c)

times :: (Num a, Ord a) => (a, a) -> (a, a) -> (a, a)
times (a,b) (c,d) = (min a b * min c d, max a b * max c d)

divide :: (Fractional a, Ord a) => (a, a) -> (a, a) -> (a, a)
divide (a,b) (c,d) = if c < 0 && d > 0 then error "divide by 0"
else (min a b / max c d, max a b / min c d)

toCenter :: Fractional a => (a, a) -> (a, a)
toCenter (a,b) = ((a+b) / 2, (b-a) / 2)

fromCenter :: Num a => (a, a) -> (a, a)
fromCenter (a,b) = (a-b, a+b)
```
3. Lautaro Pecile said
```from collection import namedtuples

class Interval (namedtuple('Interval', 'a b')):
__slots__ = ()

def __add__ (self, other):
return Interval(self.a + other.a, self.b + other.b)

def __sub__ (self, other):
return Interval(self.a - other.a, self.b - other.b)

def __mul__ (self, other):
l = [x*y for x in self for y in other]
return Interval(min(l), max(l))

def __div__ (self, other):
l = [x/float(y) for x in self for y in other]
return Interval(min(l), max(l))

@property
def toCenter(self):
a = (self.a + self.b) / 2.0
b = (self.b - self.a) / 2.0
return Interval(a, b)

@property
def fromCenter(self):
a = (self.a - self.b)
b = (self.a + self.b)
return Interval(a, b)
```
4. Graham said

My Python submission is available on codepad.org.
It requires the fractions module (new in Python 2.6), which adds support for rational number arithmetic.
My work is slightly hackish and could use better variable names, but it’ll do for now.

5. ```Interval =
Struct.new(:a, :b) do
def +(other)
self.class.new(a + other.a, b + other.b)
end

def -(other)
self.class.new(a - other.a, b - other.b)
end

def *(other)
products =
[
[ a, other.a ],
[ a, other.b ],
[ b, other.a ],
[ b, other.b ],
].map { |(c, d)| c * d }
self.class.new(products.min, products.max)
end

def /(other)
quotients =
[
[ a, other.a ],
[ a, other.b ],
[ b, other.a ],
[ b, other.b ],
].map { |(c, d)| Rational(c, d) }
self.class.new(quotients.min, quotients.max)
rescue ZeroDivisionError
raise ZeroDivisionError, "divided by #{other}", caller
end
end

# Examples:
a = Interval.new(5, 7)
b = Interval.new(0, 3)

a + b # => #<struct Interval a=5, b=10>
a - b # => #<struct Interval a=5, b=4>
a * b # => #<struct Interval a=0, b=21>
b / a # => #<struct Interval a=(0/1), b=(3/5)>

## ZeroDivisionError
a / b rescue \$! # => #<ZeroDivisionError: divided by #<struct Interval a=0, b=3>>
```
6. Axio said

I guess something like this is good enough:

(defclass interval-class ()
((min :initarg :min)
(max :initarg :max)))

(defmethod print-object ((iv interval-class) stream)
(format stream “[~a,~a]” (slot-value iv ‘min)
(slot-value iv ‘max)))

(defvar it1 (make-instance ‘interval-class :min 20 :max 30))
(defvar it2 (make-instance ‘interval-class :min 10 :max 20))
(defvar it3 (make-instance ‘interval-class :min -10 :max 20))

(defmacro it-op (name min max &optional guard)
`(defmethod ,name ((i1 interval-class)
(i2 interval-class))
(macrolet ((a () `(slot-value i1 ‘min))
(b () `(slot-value i1 ‘max))
(c () `(slot-value i2 ‘min))
(d () `(slot-value i2 ‘max)))
(unless ,guard
(make-instance ‘interval-class :min ,min :max ,max)))))

(it-op it-add (+ (a) (c)) (+ (b) (d)))
(it-op it-sub (- (a) (c)) (- (b) (d)))
(it-op it-mul
(min (* (a) (c)) (* (a) (d))
(* (b) (c)) (* (b) (d)))
(max (* (a) (c)) (* (a) (d))
(* (b) (c)) (* (b) (d))))
(it-op it-div
(min (/ (a) (c)) (/ (a) (d))
(/ (b) (c)) (/ (b) (d)))
(max (/ (a) (c)) (/ (a) (d))
(/ (b) (c)) (* (b) (d)))
(or (and (<= (a) 0) (<= 0 (b)))
(and (<= (c) 0) (<= 0 (d)))))

Next time, I’ll implement a reader to add intervals to the syntax :)

7. Axio said

With less bugs :)
Please disregard the above untested code :)

(defclass interval-class ()
((min :reader i-min :initarg :min)
(max :reader i-max :initarg :max)))

(set-macro-character #\[
(lambda (stream char)
(let* ((lst (read-delimited-list #\] stream t)))
(if (= 2 (length lst))
(make-instance ‘interval-class :min (apply #’min lst) :max (apply #’max lst))
(error “Invalid interval”)))))
(set-macro-character #\] (get-macro-character #\)))

(defmethod print-object ((iv interval-class) stream)
(format stream “[~a,~a]” (i-min iv) (i-max iv)))

(defmethod plus ((i1 interval-class) (i2 interval-class))
(let ((a (+ (i-min i1) (i-min i2)))
(b (+ (i-max i1) (i-max i2))))
(make-instance ‘interval-class :min (min a b) :max (max a b))))

(defmethod minus ((i1 interval-class) (i2 interval-class))
(let ((a (- (i-min i1) (i-max i2)))
(b (- (i-max i1) (i-min i2))))
(make-instance ‘interval-class :min (min a b) :max (max a b))))

(defmacro m-d (opname op &optional guard)
`(defmethod ,opname ((i1 interval-class) (i2 interval-class))
(unless ,guard
(let* ((b0 (loop for i in (list #’i-min #’i-max)
append (loop for j in (list #’i-min #’i-max)
collect (,op (funcall i i1) (funcall j i2)))))
(b1 (apply #’min b0))
(b2 (apply #’max b0)))
(make-instance ‘interval-class :min b1 :max b2)))))

(m-d mul *)
(m-d div / (and (<= (i-min i2) 0) (<= 0 (i-max i2))))

(defvar x [1 2])
(defvar y [4 3]) ;; reverse order \o/
(list (plus x y) (minus x y) (mul x y) (div x y))

8. slabounty said

A ruby version using the Rational class. The to_s is has a bunch of extra code to make it look “pretty” if the denominator is 1. It might be better to monkey patch this in the Rational class itself.

```class Interval

def initialize(a, b)
@a = Rational(a/1)
@b = Rational(b/1)
end

return Interval.new(@a + i.a, @b + i.b)
end

def subtract(i)
return Interval.new(@a - i.b, @b - i.a)
end

def multiply(i)
return Interval.new(
[@a*i.a, @a*i.b, @b*i.a, @b*i.b].min,
[@a*i.a, @a*i.b, @b*i.a, @b*i.b].max
)
end

def divide(i)
begin
return Interval.new(
[@a/i.a, @a/i.b, @b/i.a, @b/i.b].min,
[@a/i.a, @a/i.b, @b/i.a, @b/i.b].max
)
rescue ZeroDivisionError
return "Divide by zero"
end
end

def ends_center
return Interval.new(((@a + @b) / 2.0), ((@b - @a) / 2.0))
end

def center_ends
return Interval.new(@a - @b, @a + @b)
end

def to_s
"#{a.numerator}" + ((a.denominator !=  1) ? "/#{a.denominator}" : "" ) + ", " +
"#{b.numerator}" + ((b.denominator !=  1) ? "/#{b.denominator}" : "" )
end
end

x = Interval.new(1, 2)
y = Interval.new(3, 4)
puts "x = #{x}"
puts "y = #{y}"
puts "interval x+y = #{x.add(y)}"
puts "interval x+y = #{x.subtract(y)}"
puts "interval x*y = #{x.multiply(y)}"
puts "interval x/y = #{x.divide(y)}"
puts "ends_center x = #{x.ends_center}"
puts "center_ends x = #{Interval.new(1.5, 0.5).center_ends}"
puts "interval x/x = #{x.divide(x)}"
```
9. Khanh Nguyen said
```open System

let plus (x:float*float) (y:float*float) =
(fst x + fst y, snd x + snd y)

let minus (x:float*float) (y:float*float) =
(fst x - fst y, snd x - snd y)

let times (x:float*float) (y:float*float) =
match x, y with
| (a,b), (c,d) -> (List.min [a*c; a*d; b*c; b*d], List.max [a*c; a*d; b*c; b*d])

let div (x:float*float) (y:float*float) =
match x, y with
| (a,b), (c,d) -> (List.min [a/c; a/d; b/c; b/d], List.max [a/c; a/d; b/c; b/d])
```
10. As Chun Kin Lee pointed out on my blog, my original attempt at removing duplication from times and divide doesn’t work on ranges with negative numbers, so you actually have to get the minimum/maximum of the four possible combinations.