### Interval Plots

If we compute for an interval , we know that for . Therefore, the points of the graph of (that is, ) lie in the rectangle with lower left corner and upper right corner . This observation allows us to produce interval plots of functions. By subdividing the domain interval i into many small pieces, we can hope to close in on the values of .

Here is a simple function with interesting behavior near .

This is our first domain interval.

Here we see an ordinary plot of in this domain.

Applying to the interval gives us a rather crude idea of where the points of the graph of lie.

This definition produces a rectangle from a pair of intervals.

Here we superimpose the rectangle on the graph of .

We define a function `split` so that `split[Interval[``]]` returns the two intervals `{Interval[``], Interval[``]}`.

We split the domain into two intervals.

Now, we apply to these two intervals.

A transposition gives us the two rectangles.

Here we see that we have successfully narrowed the possible values of on the second interval.

Here we split the domain five more times to obtain intervals.

We apply to the intervals by repeating the previous computation, which was set up to be independent of the number of intervals present.

Here we see the result. In areas where varies only little, the rectangles are oblong and we could stop the subdivisions. In other areas, the rectangles are still rather tall and we should proceed with further subdivision.

The observations from our example suggest that we should subdivide the domain intervals until the length of is small enough (within the resolution of the image). Because the function being plotted could have singularities (such as 0 for ) we should also stop if the length of gets small.

This is an adaptive technique, because not all intervals will be subdivided the same number of times. Such a method avoids unnecessary work and stops when at least one of several criteria is fulfilled.

The command `IntervalPlot[`expr`,{`x`,``,``}]` implements these ideas. It is shown, together with most of the auxiliary definitions, in Listing 1. The auxiliary procedure `refine[`f`,`init`,`tol`]` contains the code for the adaptive subdivision, starting with an initial interval init. In the main loop, the nested list `new` contains the queue of intervals that have not yet been dealt with. We take the first of these intervals and compute . Note that we put the domain interval into a list. The reason for this will become clear in the next section when we use the same code for three-dimensional plots. If either or is small enough (smaller than tol, the value of the `Tolerance` option), we add the rectangle to the list of results; otherwise we split and put the resulting intervals back into the queue.

The queue is maintained as a nested list, , instead of a simple list because append and drop operations are relatively slow compared to the manipulations of a nested list. The language LISP uses the same data structure. Table 1 compares the operations of adding an element to the queue and removing one for ordinary lists and nested lists. Note that the combined operation of obtaining the first element and removing it from the list is particularly simple for nested lists (because of Mathematica's parallel assignments of list elements). We use a similar nested structure for the list of final rectangles. Because the elements of this nested list are lists themselves, we use a head `h`, different from `List`, to hold the elements together to avoid problems when flattening the list. The nested list is converted into an ordinary list with the help of `Flatten`.

The function `interval[`arg`]` deals with cases where Mathematica's interval arithmetic does not return a proper interval, either because it returns a single number (which is turned into an interval) or for some other reason, in which case we return the interval . If an interval boundary is infinite, we replace it by a huge finite value and color the rectangle red to alert us to the singularity present.

 --
Table 1. Recursive and iterative queue operations.

```
Options[IntervalPlot] = { Tolerance -> 0.01 }

length[i_] := Max[i] - Min[i]

split[i_Interval] :=
With[{m = (Min[i] + Max[i])/2},
{Interval[{Min[i], m}],Interval[{m, Max[i]}]} ]

SetAttributes[interval, Listable]

interval[i_Interval] := i
interval[i_?NumericQ] := Interval[i]
interval[_] := Interval[{-Infinity, Infinity}] (* junk *)

infcolor = RGBColor[1,0,0]; huge = 10.0^10;

rect[{xi__,yi_}] /; Min[yi] == -Infinity || Max[yi] == Infinity :=
{infcolor, rect[{xi, yi /. DirectedInfinity[d:(-1|1)] :> d huge}]}

rect[l:{_,_}] := Rectangle[Min /@ l, Max /@ l]

RectangleGraph[l:{{_,_}..}] := Graphics[{Thickness[0], rect /@ l}]

IntervalPlot[expr_, {x_Symbol, x0_, x1_}, opts___?OptionQ]:=
With[{tol = (x1-x0)Tolerance /. {opts} /. Options[IntervalPlot]},
Module[{finals},
finals = refine[Function[x, expr],
Interval /@ N[{{x0, x1}}], tol];
Show[RectangleGraph[finals],
Evaluate[FilterOptions[Graphics, opts]],
PlotRange -> {{x0, x1}, Automatic}, Axes -> True]
] ]

refine[f_, init_, tol_] :=
Module[{new = {init, {}}, finals = h[], i, j},
While[new =!= {},
{i, new} = new; (* dequeue *)
j = interval[f @@ i];
If[ length[j]<tol || Max[length /@ i] < tol,
finals = h[Append[i,j], finals]
,(*else*)
Scan[ (new = {#, new})&, split[i] ] (* queue new ones *)
]
];
List @@ Flatten[finals, Infinity, h]
]
```
Listing 1. Interval plots.

Our code is in the `Intervals.m` package which we load after restarting the kernel.

The difficult function gives a rather good result, including the singularity at zero.

For the function , the simple pole is clearly marked in red.

The function is not treated in the best possible way at zero. The reason is that the numerator and denominator of the fraction are treated as independent values, as explained in the previous section.

A function with periodic features that are smaller than the tolerance shows the effects of aliasing. (Aliasing denotes features that occur in a sampled function when the sampling rate is too low for the true shape of the function to be reconstructed; it is observed frequently in digitized images and sounds.) In an ordinary plot, we might miss essential features; here, we are alerted to the problems present.

If the tolerance is smaller than the period of the oscillation, a satisfactory picture results.

This next example is often used to justify the use of interval plots. The small singularity is clearly detected.

An ordinary plot may miss the singularity and show only a small spike.

Converted by Mathematica      October 15, 1999 [Prev Page][Next Page]