Time Space Tradeoffs in Vector Algorithms for APL Functions

Transcribed from the printed issue; not yet reviewed. Read it in the PDF of the issue, page 101

Abstract¶

When considering algorithms which will run on processors in which large vectors of values can be manipulated in time which is independent of the length of the vector, many issues of time versus space arise which are not important, or of less importance, on more conventional machines. This paper explores various alternative techniques which might be used in an implementation of APL on such a machine, and considers some of these time/space tradeoffs.

Keywords: code generation, vector machines, APL, time/space tradeoff.

Introduction¶

An observed characteristic of many APL expressions is what might be called a "bulging middle". That is, they may start out with relatively small inputs, which are then transformed (by an inner or outer product, for example), into large intermediate values. These intermediate results are then compacted, for example by reduction or compression, into relatively small outputs. The classic idiom for generating a boolean vector representing the primes less than \(n\) [Perlis & Rugaber 1979] is a good example:

2=+⌿0=(⍳n)∘.|⍳n

Here the inputs are scalar and the output is a vector of length \(n\), but in between an array of size \(n^2\) is both constructed and reduced. Many APL programmers have had the annoying experience of writing a function which works fine on small test cases, but then runs out of memory on larger data because the system could not accommodate these intermediate results.

Some APL implementations [Budd 1988] have attempted to circumvent this problem by means of a lazy evaluation approach which produces elements on demand. Using these techniques the large intermediate results are never generated all in one place, but rather are produced element by element as needed. In this manner memory requirements can be kept to a minimum. However little or no time savings are produced using these algorithms, except in the cases such as compression where portions of the computation which do not contribute to the result can be eliminated altogether.

This situation changes when one starts to consider algorithms which might be used to implement various APL functions on parallel architectures, such as vector machines. In this paper we will present examples illustrating some of the time/space tradeoffs which can occur. This is not to say there are not time/space tradeoffs which can be made in scalar implementations of APL. Indeed, for some functions, such as membership (dyadic ⍳) the alternatives are dramatic [Budd 1981]. However, for most functions, and certainly for the most commonly occurring functions, the choices are straightforward. This is not true for vector implementations.

The model we will use for describing the vector algorithms is that used by Hillis and Steele [Hillis & Steele 1986] to describe the Connection Machine [Hillis 1985]. Our notation for algorithms will be Algol-like, with the addition of an instruction:

for all \(k\) to \(n\) in parallel do \(s\) od

which will execute the statement \(s\) in parallel \(n\) times, each instance of \(s\) having access to a unique value \(k\). The total time for this computation will be assumed to be the same as the time necessary for a single execution of \(s\), assuming \(n\) does not exceed some maximum vector length (on the Connection Machine, this limit is the number of processors). For larger values of \(n\), the vector is divided into equal components each of the maximum size, and an implicit loop is made around the statement \(s\). Our notation differs from Hillis and Steele in that it explicitly provides the upper bound \(n\).

Reduction and Outer Product¶

Saal and Weis reported [Saal & Weis 1975] that reduction and outer product were among the most commonly used APL functions. These two frequently also occur in combination. A somewhat contrived example, but similar to many real situations (such as the primes expression of the last section), will illustrate a few of the possible tradeoffs that can be made. Our example expression will simply take a scalar value \(n\), produce an \(n^2\) array using outer product, then reduce this twice to once more produce a scalar result.

+/+⌿(⍳n)∘.+⍳n

The naïve and straightforward algorithm for this example would evaluate the expression from inside out, forming an entire intermediate value before considering the next outermost function. In this example this is not only the most straightforward algorithm, it is also the most time efficient. All \(n^2\) elements of the outer product can be generated in one step, as follows:

for all k to n*n in parallel do
      t[k] := (((k-1) mod n) + 1) + (((k-1) div n) + 1)
od

The innermost reduction can then be performed using a variation of Hillis’ and Steele’s array summing algorithm [Hillis & Steele 1986], which works in \(\mathrm{O}(\log_2 n)\) time:

for i := 1 to ceil(log2 n) do
    for all j to n*n in parallel do
        if ((j-1) div n) mod 2^i = 0 then
          if ((j-1) div n) + 2^(i-1) < n then
            t[j] := t[j] + t[j + (2^(i-1) * n)]
    od
od

And the outermost reduction in a similar fashion. Thus we can produce the answer using \(\mathrm{O}(\log_2 n)\) time and \(\mathrm{O}(n^2)\) space.

Supposing \(n\) is large, and we would prefer not to generate the \(n^2\) intermediate values. Using demand-driven techniques [Budd 1988], we could combine elements of the outer product with the reductions. Since the nested reduction is computing column-sums, one approach would be to use a slightly less efficient algorithm for the reduction, but combine the summing of the reduction with the generation of rows for the outer product:

for all i to n in parallel do
      t[i] := 0
od
for i := 1 to n do
      for all j to n in parallel do
            t[j] := t[j] + (i + j)
      od
od

Here we use \(\mathrm{O}(n)\) time, but only \(\mathrm{O}(n)\) space. The outermost reduction can be accomplished, as before, using an \(\mathrm{O}(\log_2 n)\) summation to produce the final scalar result.

But even this can been viewed as being space inefficient, since the final output is simply a scalar. The final possibility is to ignore vector instructions altogether, producing instead single values from the outer product on demand.

t := 0
for i := 1 to n do
    for j := 1 to n do
       t := t + (i + j)
    od
od

So we have a sequence of possible algorithms to implement the functionality of this expression:

Which is the "correct" approach? The answer, of course, is that under various circumstances any one of them might be the correct technique to use, depending upon the amount of memory available and how quickly the user wants the answer produced. One problem is that while it is easy to optimize for time (the naïve first algorithm, which simply produces intermediate results from the inside out), it is much more difficult for a system to make compromises with time and space, since these involve the interaction of adjacent functions.

It should be noted that matrix multiplication (inner product in APL terms) is in many ways similar to the reduction of an outer product. It should come as no surprise, therefore, that there are several possible choices in matrix multiplication algorithms of a vector nature, from \(\mathrm{O}(\log_2 n)\) time and \(\mathrm{O}(n^3)\) space, to \(\mathrm{O}(n)\) time and \(\mathrm{O}(n^2)\) space, for complete results. If only partial results are desired (for demand driven evaluation), there exist \(\mathrm{O}(n)\) time and \(\mathrm{O}(n)\) space algorithms, and a \(\mathrm{O}(n^2)\) time and O(1) space version.

Why would an implementation choose to use a less time efficient algorithm? In our model, if the size of a vector operation exceeds the maximum permitted size the vector is simply divided into equal portions and a loop is made around the separate pieces. If \(n^2\) is close to this maximum length, then a \(\mathrm{O}(\log_2 n)\) operation on a vector of size \(n^3\) will still take \(\mathrm{O}(n \log_2 n)\) steps. Thus one of the more space efficient algorithms would be likely preferable, despite their asympotically worse performance.

A Searching Example¶

Our second example again exhibits a time/space tradeoff, and also illustrates some of the non-obvious algorithms which can be considered when vector instructions are available that work in constant time.

The APL function "index of" (dyadic ⍳) takes two arguments, the first of which must be a vector. It returns a result of the same size and shape as the right argument, only containing, for each element in the right argument, the position of the first occurrence of the element in the vector left argument. The length of the vector left argument plus one is returned for values which do not occur. For example 27 32 32 17 ⍳ 12 17 27 32 19 yields 5 4 1 2 5.

To describe the possible algorithms which might be used to implement this function, let us suppose the left argument, call it \(a\), is a vector of length \(n\) and the right argument, \(b\), a vector of length \(m\). Furthermore assume we want to place the answer into a vector of length \(m\) called \(c\).

The naïve algorithm is simply a double loop, which iterates over the values of b, checking each against the successive values of a.

for i := 1 to n do
    c[i] := m + 1
    for j := m downto 1 do
       if a[i] = b[j] then
          c[i] := j
    od
od

Using this algorithm takes \(\mathrm{O}(n * m)\) operations to produce the \(n\) result values.

It has long been known [Bernecky 1973], that execution of this function can be speeded considerably by preprocessing the left argument, producing an array of length \(n\) containing the positions of the elements in sorted order. (This is, in fact, what is produced by the APL sorting function grade-up). For example if a is the vector 17 32 12 97, the temporary variable containing the sorted positions would be 3 1 2 4, indicating the fact that the smallest value occurs in position 3, the next largest in position 1, and so in. Using this temporary variable and three additional variables one can then use a binary search to find the position of the arguments from the right argument. Let us assume this temporary sorted value has already been generated and placed in the variable \(t\).

The algorithm to produce the \(n\) output values is therefore as follows:

for i := 1 to m do
    high := m
    low := 1
    for j := 1 to ceil(log2 n)  do
         mid := (high + low) div 2
         if a[t[mid]] < b[i] then low := mid
         else high := mid
    od
    if b[i] = a[t[mid]] then c[i] := t[mid]
    else c[i] := m + 1
od

Using this technique the time necessary to produce the results has been reduced to \(\mathrm{O}(n \log_2 n)\) for sorting the left argument, and \(\mathrm{O}(m \log_2 n)\) for performing the \(m\) binary searches, at the cost of an additional \(\mathrm{O}(n)\) space that would not otherwise be necessary (the variable \(t\)). Simply by changing the implementation of the membership algorithm in this fashion, execution time of one example expression was observed to drop from 34.2 seconds down to 2.0 seconds [Budd 1981].

When one considers the possibility of implementing the membership function on a vector machine, it is tempting to assume that the preprocessing of the left argument is no longer necessary. Instead, in one vector step it is possible to discover whether any position matches a given element in the left argument; for example as in the following (which assumes no elements in the left argument are repeated):

for all i to n in parallel do
      c[i] := m + 1
od
for i := 1 to n do
      for all j to m in parallel do
            if a[i] = b[j] then c[i] := j
      od
od

The answer is produced in \(\mathrm{O}(n)\) steps and requires no additional space. In working with vector algorithms, however, one quickly learns to be suspicious of linear time algorithms. And indeed, by once more preprocessing the left argument (which now can be done in \(\mathrm{O}(\log_2 m)\) vector steps [Hillis & Steele 1986]), we can perform all \(n\) binary searches in parallel, producing the answer in \(\mathrm{O}(\log_2n)\) time. (But with a cost of four size \(n\) temporary variables.)

for all i to n in parallel do
      high[i] := m
      low[i] := 1
od
for j := 1 to ceil(log2 n)  do
      for all i to n in parallel do
            mid[i] := (high[i] + low[i]) div 2
            if a[t[mid[i]]] < b[i] then low[i] := mid[i]
            else high[i] := mid[i]
      od
od
for all i to n in parallel do
      if b[i] = a[t[mid[i]]] then c[i] := t[mid[i]]
      else c[i] := m + 1
od

Conclusions¶

In moving from an environment in which values are produced one at a time on a scalar processor to one in which entire vectors can be manipulated in constant time, not only do many of our conventional views and opinions about algorithms become invalid, but choices and alternatives arise which are not present, or not as critical, in the scalar world. In particular there are time versus space tradeoffs which become much more pronounced in the vector environment than they were in the scalar case. Of course, from a theoretical point of view we should not be surprised at all that we can trade polynomial time and constant space (as in the reduction example) for logarithmic time and polynomial space [Parberry 1986]. The question is not whether such tradeoffs are possible, but under what circumstances the choice should be made.

In this paper we have presented a number of alternative schemes for implementing several operations in the language APL. Further work will be necessary to develop heuristics to help the language implementor in selecting the most appropriate technique.

Timothy A. Budd, Department of Computer Science, Oregon State University
September 20 1994

References¶

[Bernecky 1973]

Bernecky, B., "Speeding Up Dyadic Iota and Dyadic Epsilon", In Proceeding of the APL Congress 73, (Copenhagen, Denmark, Aug 22-24). Elsevier North-Holland, New York, 1973, pp. 479-482.

[Budd 1981]

Budd, Timothy A., "An APL Compiler", Technical Report 81-17, University of Arizona, October 1981.

[Budd 1988]

Budd, Timothy A., An APL Compiler, Springer-Verlag, 1988.

[Hillis 1985]

Hillis, W.D., The Connection Machine, MIT Press, Cambridge, Mass. 1985.

[Hillis & Steele 1986]

Hillis, W. D. and Steele, G. L., "Data Parallel Algorithms", Communications of the ACM, Vol 29(12):1170-1183 (December 1986).

[Parberry 1986]

Parberry, I., "Parallel Speedup of Sequential Machines: A Defense of the Parallel Computation Thesis", Sigact News, Vol 18(1), 1986, pp. 54-67.

[Perlis & Rugaber 1979]

Perlis, Alan J. and Rugaber, Spencer, "Programming with Idioms in APL", APL Quote Quad, Vol 9(4):232-235 (June 1979).

[Saal & Weis 1975]

Saal, H. and Weis, Z., "Some Properties of APL programs", APL ’75 Conference Proceedings, (Pisa, Italy, June 11-13). ACM, New York, 1975, pp. 292-297.

This article was originally written in 1994. We thank Professor Budd for permission to print it in Vector - Ed.