interval arithmetic & assume database



I wasn't aware that Mathematica had such interval capabilities.  Interesting!  Is there a paper somewhere that describes this system?

Re complex "regions":

The Apple Mac has had a scheme for handling "arbitrary" 2D regions since almost day one.  This sort of capability is required for many window & graphics operations; e.g., you need a "mask" for each arbitrary-shaped icon.  Mac regions are completely arbitrary portions of the plane, including those with holes and infinite (?) regions.  I think that Mac regions are bitmap-based, so there is an inherent resolution limitation.

Another scheme for handling arbitrary 2D regions is based on quadtrees, so that portions of regions requiring greater resolution can be modeled more precisely.

One could triangulate the plane, but I'm not sure what benefits this scheme would have over using rectangular regions.

Purely convex regions can be modelled by a multiplicity of linear boundary lines/line-segments, but this does not produce a Boolean algebra.

One problem with the plane is that it is no longer so easy to represent open & closed edges, corresponding to (half-)open and (half-)closed intervals.  For example, trying to represent exactly which points along the circle defined by the radius of convergence of a power series might be difficult.

It would appear that any more sophisticated scheme for boundary curves -- e.g., circles, ellipses, conic sections, cubic splines, etc. -- will have severe difficulties dealing with unions, intersections, and complements.

One might take inspiration from my "gradual precision" intervals in order to handle 2D regions with curved boundaries.  Graphics programs utilize (typically rectangular) "bounding boxes" to make quick one-way decisions about region intersections, and if the quick decision returns "not sure", then one has to go back & evaluate the bounding curve to greater precision -- e.g., by doing another several levels of quadtree recursion.

This "gradual precision" method for disentangling different regions can be arbitrarily difficult -- e.g., the separation of the complex plane into 4 distinct regions by the complex Newton method for finding the roots of z^3-1, where the four regions are: the points that eventually converge to 1, the points that eventually converge to %e^(2*%pi/3), the points that eventually converge to %e^(-2*%pi/3), and the points that don't converge at all.  These regions are classical fractal patterns.

Inspiration can also be taken from so-called "vector graphics" formats (as opposed to "bitmap graphics" formats): SVG, OpenGL, etc.

One can, of course, create a pseudo-Boolean algebra by simply composing the mathematical expressions describing certain regions together with AND's, OR's and NOT's, but all we are left with is a procedure for deciding if a particular point is inside the region or not.  We cannot, in general, reduce these expressions to a canonic form to decide equivalence, nor can we decide if two such regions are disjoint.

At 03:40 PM 2/4/2013, David Stoutemyer wrote:
>I like all of these suggestions, and I am glad to see that they include both open and closed endpoints.
> 
>The suggestion for representation as multi-intervals, with disjoint ordered unions of intervals is realized in Mathematica, and most of its functions seem to know how to handle such multi-intervals, though not as thoroughly  their adaptive significance arithmetic.
> 
>An additional area to consider is complex multi-intervals, even though they might entail some of the problems mentioned for multivariate real multi-intervals. I think that some work has been done (by W. Kahan?) using circular or ellipsoidal regions to reduce growth, but in any event people would surely at least want multi-intervals that are unions of rectangular intervals. One representation would be a real multi-interval + %i times another real multi-interval. Another representation would be as a lexically ordered union of rectangles each specifed by its lower left and upper right corners. The implementation would be simpler if simplification produces only one of these alternatives. Otherwise numerous complicated rules must be implemented to simplify combinations of the two forms well.
> 
>It would also be nice if users could use rewrite rules to add multi-interval handling to their own functions or those of others  without doing surgery on the built-in interval implementation or other peoples function implementations.
> 
>I would also like to see the infectiousness of Floats and Big Floats be on a local basis of each real arithmetic operation rather than degrading all numbers in a multi-interval whenever one of them is a Float or a Big Float.
> 
>-- best regards, david stoutemyer
>
>On Mon, Feb 4, 2013 at 12:26 PM, Henry Baker <hbaker1 at pipeline.com> wrote:
>After all of the discussion over the past several months about interval arithmetic & the assume database, I got to thinking about what would be an appropriate fragment of a decidable logical theory that one could easily implement in Maxima, and also would be efficient enough for general usage.
>
>Below is my suggestion for such a "decidable interval theory".
>
>Rather than presenting axioms and rules of inference, I will go directly to a model.
>
>We have _representable real numbers_ ("RRN"), including extensions with +oo (positive infinity) and -oo (negative infinity), together with three comparison predicates: =, <, <=.  These predicates have the usual trichotomy property aka "total/linear ordering".
>
>Representable real numbers include the rational numbers, positive n-th roots of positive representable reals, certain transcendental constants (%pi, %e), etc.
>
>We start with "atomic" intervals: for any x<y in RRN, we have 4 distinct atomic intervals: [x,y], [x,y}, (x,y], (x,y), with their usual meanings.  We let the symbol "{" stand for either "[" or "(", and the symbol "}" stand for either "}" or ")"; thus {x,y} is any of the 4 intervals [x,y], [x,y), (x,y], (x,y).
>
>For any x in RRN, we have the singleton set [x].  By an abuse of notation, we can also write [x,x]=[x].  The notation {x,x} can mean only the singleton [x].
>
>Given two atomic intervals {w,x}, {y,z}, we can determine if they are _disjoint_.  W/o loss of generality, we assume w<=y.
>
>{w,x} is disjoint from {y,z} if
>1) x<y or
>2) x=y and {w,x}={w,x) and {y,z}=(y,z}
>
>Given two non-disjoint atomic intervals {w,x} and {y,z}, we can compute their _intersection_ as sets; this intersection will necessarily also be an atomic interval:
>
>1) If w<y and z<x, then {w,x} encloses {y,z}, so the intersection is {y,z}.
>2) If y<w and x<z, then {y,z} encloses {w,x}, so the intersection is {w,x}.
>...
>
>We can now introduce _finite unions of mutually disjoint atomic intervals_.  The empty union consists of no atomic intervals at all, and is hence the empty set.  We extend the notation such that [x,x) = (x,x] = (x,x) = empty set.  It is convenient if the finite union is ordered by increasing order of the first element (implies also ordered by the second element).
>
>Given an atomic interval {x,y}, we can compute the _complement_ of the interval as follows:
>
>1) {x,y} = [x,y], then the complement is [-oo,x) union (y,oo];
>2) {x,y} = (x,y], then the complement is [-oo,x] union (y,oo];
>3) {x,y} = [x,y), then the complement is [-oo,x) union [y,oo];
>4) {x,y} = (x,y), then the complement is [-oo,x] union [y,oo].
>
>If we have two ordered finite unions of mutually disjoint atomic intervals, it is obvious how to _merge_ the ordered finite unions to produce the union of all of the intervals, coalescing atomic intervals where possible.  This process has complexity linear in the total number of atomic intervals.
>
>Also, if we have an ordered finite union of mutually disjoint atomic intervals, it is also obvious how to _complement_ this union to produce an ordered finite union of disjoint atomic intervals which represents the complement of this set.  This process has complexity linear in the total number of atomic intervals.
>
>We now have a Boolean algebra whose elements are finite unions of mutually disjoint atomic intervals.
>
>A note on performing ordering comparisons of irrational & transcendental numbers.  Since we allow the representation of certain irrational & transcendental numbers, it will be necessary to compare arbitrarily complex expressions for ordering purposes.  Since we have no way of knowing in advance how many digits of precision will be required to determine the ordering of these arbitrarily complex expressions, we will need to represent these numbers by an object containing the algebraic expression itself, together with a pair of rational numbers (possibly including oo and -oo) which have the current tightest bounds on the irrational/transcendental number.  If/when it becomes necessary to _refine_ the precision of these bounds in order to perform an accurate comparison, the algebraic expression will be evaluated to higher precision to compute tighter bounds.
>
>We are finally in a position to utilize this Boolean algebra in Maxima.
>
>Our "assume database" will contain for each variable in the database an association with an element of the Boolean algebra.  Each "real" variable will initially be associated with the atomic interval (-oo,oo).  We denote the Boolean algebra element associated with x in assume database D by B(x,D).
>
>If Maxima encounters an expression of the form (if x<c then E1 else E2) in an environment where the current assume database is D, then the expression E1 will be simplified assuming that x is associated with intersect(B(x,D),(-oo,c)), and the expression E2 will be simplified assuming that x is associated with intersect(B(x,D),[c,oo)).
>
>Equality is handled as x=c <=> x<=c and c<=x; i.e., x is associated with [c].
>
>If Maxima encounters an expression of the form (if x<c1 and y<c2 then E1 else E2), then this is equivalent to
>
>(if x<c1 and y<c2 then E1
> elsif x<c1 and not y<c2 then E2
> elsif not x<c1 and y<c2 then E2
> elsif not x<c1 and not y<c2 then E2)
>
>We will also need a rule whereby equal arms can be coalesced back together: (if ... then E2 else E2) => E2
>
>Note that our "decidable interval theory" cannot handle relations between two different variables, e.g., (if x<y then E1 else E2).  (It is possible to develop a decidable theory of such n-dimensional "regions", but the computational complexity of this theory explodes; even a two-dimensional theory necessary to bound regions in the complex plane could quickly become unwieldy.)
>
>However, we can _approximate_ or _bound_ the results from expressions such as x<y, but these bounds will most likely not be tight; e.g., if x is associated with (-oo,oo) and y is associated with (-oo,oo), then even if we know that x<y, this does not further constrain either x or y _by itself_.
>
>It is easy to see how the operations +,-,*,/ can be extended to operate on elements of our Boolean algebra.  Unlike traditional "interval arithmetic", our finite unions of atomic intervals is closed under division.
>
>Other functions can also be extended to operate on finite unions of atomic intervals, although the complexity goes up substantially if these functions are not strictly monotonically increasing (or at least non-decreasing).
>
>Once we have our language for discussing atomic intervals and finite unions of atomic intervals, we are in a position to better characterize _partial functions_, i.e., functions not defined over the entire real interval (-oo,oo).
>
>Thus, we can define a restricted function sin(x) which operates only on the interval -%pi/2 to %pi/2, and hence is monotonicly increasing.
>
>We can also better characterize the range of functions such as cosh(x) which is [1,oo).
>
>We can now extend Maxima with two new functions on expressions: Domain(expression) and Range(expression).  These functions are merely _approximations_ of the domain and range of the expression, but for some simpler expressions, they can be accurate.
>
>Thus, Domain(1/x) = (-oo,0) union (0,oo) = Range(1/x).
>
>Domain(sin(x)) = (-oo,oo); Range(sin(x)) = [-1,1].
>
>Actually, Domain() and Range() accept a second argument, which is a current assume database, which is used to bound the various variables within the expression.
>
>Domain(if x<c then E1(x) else E2(x),D) = Domain(E1(x),B(x,D) intersect (-oo,c)) intersect Domain(E2(x),B(x,D) intersect [c,oo))
>
>Range(if x<c then E1(x) else E2(x),D) = Range(E1(x),B(x,D) intersect (-oo,c)) union Range(E2(x),B(x,D) intersect [c,oo)).