# The cube connected cycles: a versatile network for parallel computation

Franco P. Preparata, Jean Vuillemin

<figure><img src="image_1.png" /></figure>

## To cite this version:

Franco P. Preparata, Jean Vuillemin. The cube connected cycles: a versatile network for parallel computation. [Research Report] IRIA-RR-356, IRIA. 1979, pp.35. hal-04716599

HAL Id: hal-04716599

https://inria.hal.science/hal-04716599v1

Submitted on 1 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

laboria

laboratoire de recherche
en informatique
et automatique

<figure><img src="image_4.png" /></figure>

THE CUBE-CONNECTED CYCLES :
A VERSATILE NETWORK
FOR PARALLEL COMPUTATION

Franco P. PREPARATA
Jean VUILLEMIN

Rapport de Recherche N° 356

Juin 1979

Institut de Recherche
d'Informatique
et d'Automatique

Domaine de Voluceau
Rocquencourt
B. P. 105 78150 - Le Chesnay
France
Tél.: 954 90 20
# THE CUBE-CONNECTED-CYCLES : A VERSATILE NETWORK FOR PARALLEL COMPUTATION

Franco P. Preparata*, Jean Vuillemin**

## Résumé :

Nous décrivons une structure d'ordinateur, le CCC, réalisable grâce à la technologie des circuits intégrés à grande échelle. Cette structure est programmée en détail pour traiter efficacement de divers problèmes, en particulier de la transformée de Fourier rapide, du tri, et des algorithmes dérivés. Le CCC peut aussi être utilisé en temps que structure générale d'interconnection de microprocesseurs en vue de traitements parallèles.

## Abstract :

We introduce a network structure, the CCC, complying with the present technological constraints of VLSI design. We describe in detail how to program the network for efficiently solving a class of problems which includes Fast Fourier Transform, Sorting, and derived algorithms. The CCC structure can also be used as a general purpose parallel processor.

<footnote>* On leave from the University of Illinois, Urbana, ILL 61801 (U.S.A.)</footnote>

<footnote>**Université Paris-Sud, Laboratoire de Recherche en Informatique, Bat. 490 - 91405 Orsay (France)</footnote>
# 1 - INTRODUCTION

The great technological progress embodied by VLSI has made possible to conceive large systems of processing elements cooperating in the execution of parallel algorithms. This has motivated considerable research interest in parallel computation. Unfortunately, here the situation is very different from that of serial computation, where the RAM machine [1] represents a universally accepted model. To bypass the difficult problem of developing a realistic model -thereby accepting all of its constraints- several workers [2-5] have preferred to design and analyze parallel algorithms with reference to a system where each pair of processor is connected (or an equivalent system). Although not without merit, because it aims at uncovering the inherent data-dependence of given problems, this approach ignores the technological constraints of VLSI, particularly as regards the communication among the processing elements [6]. At the opposite, other researchers [7-11] suggest that processor interconnection should be limited to planar links between topologically neighboring cells (arrays or meshes). Such designs are certainly well suited for current VLSI technology, and they have cleverly been used in implementing matrices or graph problems [9-12] for example. This type of connection is, however, not suited for efficiently implementing various fundamental problems, such as sorting and convolution. Indeed, good algorithms for solving these problems intrinsically require data movement between processors which are topologically far apart; for example, sorting on a $n \times n$ array processor such as ILLIAC IV requires time $\Omega(\sqrt{n})$ [8].

The purpose of the paper is to propose and analyze a new interconnection of processors, called the *cube-connected-cycles*, which is remarkably suited for implementing efficient algorithms such as Fast Fourier Transform (FFT), Sorting, etc.... The geometric structure underlying the interconnections is that of the $k$-dimensional cube. This structure which has already been studied in relation to parallel computation [12], is not readily usable for VLSI design, since each of the $2^k$ processors is connected to $k$ other processors.

By combining parallelism and pipelining we are able to achieve the following results :

(1) The number of connections per processor is reduced to 3 ;

(2) Processing time is not significantly increased with respect to that achievable on the k-cube structure ;
(3) Programs for the individual modules are obtained in a systematic way from a standard recursive description of the global algorithms;

(4) The overall structure complies with the basic requirements of VLSI implementation: modularity, ease of layout, simplicity of communication among the processing elements, simplicity in timing and control of the entire system [14]. We also propose a wire layout of the CCC, which can be physically realized with two layers of silicon. One layer is devoted to vertical wires, one to horizontal ones.

(5) Finally we are able, without resorting to any drastic departure from classical algol like languages, to provide fully accurate and hopefully easily understandable descriptions of our parallel programs. This is a favorable sign that VLSI design may possibly be endowed with suitable high level programming languages.

## 2 - DESCRIPTION OF THE SCHEME

Our parallel computing system, the *cube-connected-cycles (CCC)* is a network of identical processors. Each of these processors called, a *module*, contains an operand register T, a few memory locations, and possesses basic arithmetical and logical capabilities. It is controlled by a stored program. For the time being, we make the hypothesis of *unlimited parallelism*, that is, the number of modules is tailored to the problem size; under this hypothesis, the one or two memories mentioned earlier suffice. Subsequently, under the hypothesis of *limited parallelism*, we shall endow each module with a small private RAM. In either case, each module is somewhat simpler than a current microprocessor but not basically different from it.

A module has 3 interconnection ports. Each interconnection line linking two modules can be used for the bidirectional transmission of one operand, and it is irrelevant whether operand transmission is serial or parallel. For correctly executing the algorithms described in the following sections, it is indifferent to synchronize the entire system through a central clock, which defines time units for all modules, or to let synchronization problems be settled at the level of each communication line, thus achieving an asynchronous system. In order to describe the interconnections, we assume for simplicity that n, the number of modules, is a power of two, i.e. $n=2^k$; moreover, assume that k is of the form $k=r+2^r$; the modifications resulting when k is arbitrary are straightforward (in the latter case, r is the smallest integer for which $r+2^r \ge k$). Thus, for example, n=64 for r=2, n=2048 for r=3, and n=1053576 for r=4. Each module
<figure><img src="image_5.png" /><figcaption>- Figure 1 -</figcaption></figure>

has a k-bit address m which in turn is expressed as a pair ($\ell, p$) of integers represented with (k-r) and r bits respectively, such that $\ell.2^r+p=m$.

As mentioned earlier, each module has three ports : F,B, and E (mnemonic for forward, backward, external), whose connection is entirely determined by the module address ($\ell,p$), that is :

F($\ell,p$) is connected to B($\ell,(p+1)\mod 2^r$)

B($\ell,p$) is connected to F($\ell,(p-1)\mod 2^r$)

E($\ell,p$) is connected to E($\ell+\epsilon 2^p,p$)

where $\epsilon=1-\text{BIT}_p(\ell)$. Here, $\text{BIT}_p(\ell)$ is the coefficient of $2^p$ in the binary expansion of $\ell$. The interconnection scheme is displayed in figure 1. In words, the modules are grouped into $2^{k-r}$ *loops*, each loop consisting of $2^r$ modules, cyclically connected by the F-B lines. The loops are in turn interconnected as a (k-r)-cube ; if $\langle x_0,x_1,\dots,x_{k-r-1} \rangle$ are the dimensions of the (k-r)-cube, all edges along dimension $x_i$, called collectively *sheaf i*, link modules whose addresses are (.,i).
The total number of interconnection links is $3.2^{k-1} = \frac{3}{2} \cdot n$. In addition, each module communicates with the outside world by means of an extra port which will not enter our further discussions; for example we may consider the network just outlined as the operating unit of a system which processes vectors of n operands. Each module contains one operand, denoted as $T[m]$, where m is the module address.

It is interesting to examine the just described CCC within the framework of the "VLSI model of computation" recently proposed [14,18]. In this model, each wire has unit width on the silicon chip and transmits a unit of information in a unit of time; information is taken from or delivered to special areas on the chip, called nexuses, each associated with a module. Within this model, which takes realistic account of the placement of modules and interconnection, C.D. Thompson has studied the implementation of the Fast-Fourier Transform [18] and has elucidated significant relationships between chip surface A, processing time T, and a graph-theoretic quantity called the minimum bisection width $\omega$. The latter is defined for a graph $G=(V,E)$ as the minimal integer such that $\omega=|\{(u,v) \in E: u \in V_1, v \in V_2\}|$, where $\{V_1, V_2\}$ is a partition of V with $|V_1| \le |V_2| \le |V_1|+1$. Thompson has shown that, for the n-point FFT, $T \ge n/2\omega$ and $A \ge \omega^2/4$, whence $AT^2 \ge n^2/16$. With respect to the CCC interconnection pattern, we are able to show that $T=O(\log n)$ (see section 4.3) and $A=O((n/\log n)^2)$, thus achieving the bound exactly; in other words, the CCC is optimal for the FFT in the VLSI model. (The best known previous realization [19,20] attained $T=O(\log n)$ and $A=O(n^2/\sqrt{\log n})$.)

To show that $A=O((n/\log n)^2)$, we consider a layout which uses two sheaves of evenly spaced wires, horizontal and vertical. Specifically, in the $(x,y)$ plane, horizontal lines $y=i$ ($n \ge 0$) are used for cube connections, while vertical lines $x=j$ ($j \ge 0$) are used for loop connections; of the latter, lines $x=2h$ (for $h \ge 0$) are used for the return connections of the loops. Ignoring momentarily the wire U-turn necessary to form the loops (which can be achieved by extending the chip by 2 units in the y-direction), we note that a $(2^2-1) \times (2.2^2)$ chip realizes the CCC-interconnection for $2.2^2$ modules (see figure 2a); assuming inductively that the CCC for $(j-1)2^{j-1}$ modules is realizable on a $(2^{j-1}-1) \times (2.2^{j-1})$ chip, we conclude that $j \cdot 2^j$ module CCC can be placed on a $(2^j-1) \times (2.2^j)$ chip (see figure 2b). Thus, letting $n=j2^j$ we obtain $j \approx \log_2(\frac{n}{\log_2 n})$, whence the required chip size is about $(\frac{n}{\log_2 n} + 2) \times (2 \cdot \frac{n}{\log n}) = O((n/\log n)^2)$. A slightly more complicated construction (suggested in figure 2c) cuts the chip area by 25%.
<figure><img src="image_6.png" /><figcaption>Figure 2a : Interconnection of 2.2² modules.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 2b : Interconnection of 4.2⁴ modules.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Figure 2c : A more economical design.</figcaption></figure>
As for the dual ASCEND :

```
proc ASCEND(S,J)
```

```
if J≥0 then
```

```
pardo ASCEND(S,J-1), DESCEND(S+2^J, J-1) odpar;
```

```
foreach M: S≤M<S+2^J pardo OPER(M,J; T[M], T[M+2^J]) odpar
```

```
fi
```

```
corp ASCEND.
```

We expect our algol-like description of parallel programs to be self evident;
sequential iteration is controlled by a loop : **for** I=A **step** S **until** I=C, while
parallel operations are controlled by **foreach** I:p(I) **pardo** OP odpar.

As a general rule, we use ";" to denote sequencing of operations while instruc-
tions separated by "," are to be executed in parallel.

To clarify the duality between ASCEND and DESCEND, consider an integer

m = Σ_{0≤i<k} BIT_i(m) 2^i in binary, and define m̃ as m̃ = Σ_{0≤i<k} BIT_i(m) 2^{k-i-1},

that is, the binary representation of m̃ is obtained by **reflection** from that of m.
The so-called odd-even distribution (OED) [15,16] is the (involutory) permutation
computed by (see figure 2) : **for each** M:0≤M<2^k **pardo** T[M] ↔ T[M] **odpar**.

<figure><img src="image_9.png" /><figcaption>**Figure 2 :** An eight element sequence and its OED distribution</figcaption></figure>

By first applying the OED distribution to its inputs, and ASCEND algorithm
can be transformed into an equivalent DESCEND algorithm (figure 3) whose basic
operation OPER is related to the original OPER by :

$$ \tilde{OPER}(M,J: U,V) = OPER(\tilde{M},K-1-J; U,V). $$

<figure><img src="image_10.png" /><figcaption>**Figure 3 :** An algorithm and its equivalent dual.</figcaption></figure>
### 3.2) Implementation on the CCC

In order to implement ASCEND or DESCEND on the CCC, we prune the k-cube so as to use only existing connections. The first stage consists in removing the sheaves corresponding to dimensions 0,1,...,r-1, and using instead the loop connections F and B, as described in section 2. Our original (non recursive) ASCEND program is thus transformed into :

```
foreach L:0≤L<2^K-R pardo ALOOPOPER(L) odpar ;
```

```
for J←R step 1 until K-1
```

```
do foreach M:0≤M<2^K pardo if BIT_j(M)=0 then OPER(M,J; T[M], T[M+2^j]) fi odpar.
```

Here, procedure `ALOOPOPER(L)`, processes the data within loop $l$ to compute the desired result in $O(2^K)$ parallel steps, as we show later on in section 3.3. Note that the running time is still $O(k-r)+O(2^r) = O(\log n)$.

The second transformation consists in removing, for all $j=0,...,k-r-1$, the k-cube links pertaining to sheaf $(r+j)$, except those existing between modules whose addresses are of the form $(.,j)$ : the resulting interconnection is then exactly the one of CCC, as described in section 2.

The computation corresponding to the `for` loop of the above algorithm can no longer be performed in one parallel step. Using repeated cyclic shifts within loops, however, each operand in the loop can be successively brought to reside for one time unit in module $(.,j)$, where $OPER(.,J,..)$ can then be executed. Although the execution of $OPER(.,J,..)$ for all operands in a loop now requires $2^r$ time units, this computation can be *pipelined* (overlapped) with the analogous operations $OPER(.,J,..)$ for $r≤j<k$. To achieve pipelining thus requires a new function `FSHIFT(L)$$, which performs a cyclic forward shift of the operands in loop L, that is :

```
foreach J:0<J≤2^K pardo T[2^r.L+(J mod 2^K)]←T[2^r.L+J-1] odpar.
```
Our final version of `ASCEND` is thus :

`proc ASCEND`

`foreach L:0≤L<2<sup>K-R</sup> pardo ALOOPOPER(L) odpar:`

`comment` end of loop section, sheaves 0,1,..,r-1;

`for I ← 2<sup>R</sup> + 1 step 1 until I = 2<sup>R</sup>`

`do foreach L:0≤L<2<sup>K-R</sup>`

`pardo foreach P:max(I,0)≤P<min(2<sup>R</sup>,2<sup>R</sup>+I)`

`pardo if BIT<sub>p</sub>(L)=0 then OPER(A,B; U,V)`

`where A=L.2<sup>R</sup>+(P-I+1)mod 2<sup>R</sup>),`

B=P+R,

U=T[L.2<sup>R</sup>+P],

V=T[(L+2<sup>P</sup>)2<sup>R</sup>+P].

fi

`odpar;`

`FSHIFT(L)`

`odpar`

od

`corp ASCEND.`

The inner operation of the `for` loop is executed in two time units : one for `OPER`, then one for `FSHIFT`. The total running time is thus $4.2^R$ plus the time for executing `ALOOPOPER`. If we can ensure that `ALOOPOPER` can be processed in time $O(2^R)$, linear in the loop-size, the entire procedure is thus executed on the `CCC` without significant degradation of performance, in time $O(\log n)$.

A dual discussion can be presented for `DESCEND`, leading to the program :

`proc DESCEND`

`for I ← 2<sup>R</sup> - 1 step -1 until I = -2<sup>R</sup>`

`do foreach L:0≤L<2<sup>K-R</sup>`

`pardo foreach P:max(I,0)≤P<min(2<sup>R</sup>,2<sup>R</sup>+I)`

`pardo if BIT<sub>p</sub>(L)=0 then OPER(A,B;U,V)`

`where A=L.2<sup>R</sup>+(P+I-1)mod 2<sup>R</sup>),`

B=P+R,

U=T[L.2<sup>R</sup>+P]

V=T[(L+2<sup>P</sup>)2<sup>R</sup>+P].

fi

`odpar`

`BSHIFT(L) comment backwards shift of loop L;`

`odpar`

`od; comment end of treatment on sheaves k-1, k-2, ..., r;`

`foreach L:0≤L<2<sup>K-R</sup> pardo DLOOPOPER(L) odpar`

`corp DESCEND`
Figure 4 provides a pictorial view of the ASCEND procedure, both on the 3-cube and on the CCC.

<figure><img src="image_11.png" /><figcaption>Figure 4.</figcaption></figure>

### 3.3) Computation within the loops (the LOOPOPER procedures)

The next question to be answered is the implementation of LOOPOPER (L) so that it runs in time linear in the loop length. Obviously, we are constrained to using only the F and B loop links existing in the CCC.

Our objective is to emulate on the loop of length $2^r$, that is, by using only the loop links, the operation OPER as it would be executed on hypothetical r-cube sheaves. Since OPER may take place in the loop only between adjacent modules, the key permutation which realizes the desired adjacencies (with reference to ALOOPOPER) is the so-called *perfect unshuffle* [6,17]. Specifically `UNSHUFFLE(L,I)` transforms `T[L],T[L+1],...,T[L+2^{I+1}-1]` into `T[L],T[L+2],...,T[L+2^{I+1}-2],...,T[L+1],T[L+3],...,T[L+2^{I+1}-1]`, and is realized as follows:
```
proc UNSHUFFLE(L, I)
    for B + 2^I step -1 until B=2
        do foreach M:M=L.2^R+(2S+1)2^I+P with 0≤S<2^(R-I-1), -B<P<B, P mod 2 =B mod 2
            pardo T[M-1] ↔ T[M] odpar
        do
            corp UNSHUFFLE.
```

The action of UNSHUFFLE(0,I) may also be described as follows. Suppose that $T[M]=M$, for $0≤M<2^I$, and consider the contents of the loop modules as a $2^I×I$ binary address array; initially, this array contains the binary representations of the integers $[0,2^I[$ in their natural order. The array resulting after the application of UNSHUFFLE(0,I) can be obtained from the initial array by cyclically shifting its columns one position to the left (figure 5). If after UNSHUFFLE(0,I) we perform UNSHUFFLE(0,I-1), the resulting address array is obtained from the previous one by cyclically permuting the leftmost (I-1) columns, leaving the rightmost one unaffected. Thus we see that by means of unshuffle operations we can transform the initial array into one whose columns are in reverse order, i.e., we can realize OED by means of successive UNSHUFFLES:

```
proc OED(L)
    for I ← R-1 step -1 until I=1 do UNSHUFFLE(L,I) od
corp OED.
```

<table><caption>**Figure 5 :** Illustration of the action of UNSHUFFLE on the array $T[M]=M$</caption><tbody><tr><td>0 0 0</td><td>0 0 0</td><td>0 0 0</td></tr><tr><td>0 0 1</td><td>0 1 0</td><td>1 0 0</td></tr><tr><td>0 1 0</td><td>1 0 0</td><td>0 1 0</td></tr><tr><td>0 1 1</td><td>UNSHUFFLE(.,3)</td><td>1 1 0</td><td>UNSHUFFLE(.,2)</td><td>1 1 0</td></tr><tr><td>1 0 0</td><td>—&rarr;</td><td>0 0 1</td><td>—&rarr;</td><td>0 0 1</td></tr><tr><td>1 0 1</td><td></td><td>0 1 1</td><td></td><td>1 0 1</td></tr><tr><td>1 1 0</td><td></td><td>1 0 1</td><td></td><td>0 1 1</td></tr><tr><td>1 1 1</td><td></td><td>1 1 1</td><td></td><td>1 1 1</td></tr></tbody></table>

We are now ready to elucidate the general format for ALOOPOPER, which consists of a sequence of unshuffle-operation pairs, each emulating a sheaf operation. Since, however, at the end of this phase, the results are in the order induced by OED, an additional (involuntary) OED must be performed to restore the correct arrangement (see figure 6).
proc ALOOPOPER(I)
    for J<0 step 1 until J=R-1
        do foreach Q: 0≤Q<2^R, BIT₀(Q)=0
            pardo OPER(A,J; T[L·2^R+Q], T[L·2^R+Q+1])
            where A=L·2^R+(Q mod 2^J) + (Q̃ mod 2^(R-J))·2^J.
            odpar;
        UNSHUFFLE(L,R-J-1)
    od
end

OED(L)

corp ALOOPOPER

<table><caption>**Figure 6 :** A schematic presentation of ALOOPOPER for r=3.</caption><tr><td>0</td><td>1</td><td>2</td><td>3</td><td>4</td><td>5</td><td>6</td><td>7</td><td>OPER</td></tr><tr><td>0'</td><td>1'</td><td>2'</td><td>3'</td><td>4'</td><td>5'</td><td>6'</td><td>7'</td><td>UNSHUFFLE(.,3)</td></tr><tr><td>0'</td><td>2'</td><td>4'</td><td>6'</td><td>1'</td><td>3'</td><td>5'</td><td>7'</td><td>OPER</td></tr><tr><td>0"</td><td>2"</td><td>4"</td><td>6"</td><td>1"</td><td>3"</td><td>5"</td><td>7"</td><td>UNSHUFFLE(.,2)</td></tr><tr><td>0"</td><td>4"</td><td>2"</td><td>6"</td><td>1"</td><td>5"</td><td>3"</td><td>7"</td><td>OPER</td></tr><tr><td>0'''</td><td>4'''</td><td>2'''</td><td>6'''</td><td>1'''</td><td>5'''</td><td>3'''</td><td>7'''</td><td>UNSHUFFLE(.,3)</td></tr><tr><td>0""</td><td>2""</td><td>1""</td><td>3""</td><td>4""</td><td>5""</td><td>6""</td><td>7""</td><td>UNSHUFFLE(.,2)</td></tr><tr><td>0""</td><td>1""</td><td>2""</td><td>3""</td><td>4""</td><td>5""</td><td>6""</td><td>7""</td><td>OED</td></tr></table>

The execution of ALOOPOPER comprises r OPER steps, plus twice the time for OED, for a total $O(2^r)$.

The seemingly complicated expression of the parameter A in the description of A LOOPOPER(L) deserves some explanation. We recall that, for given J, A denotes the original address of the operand brought in module (L,P) by the sequence UNSHUFFLE(.,R), UNSHUFFLE(.,R-1),...,UNSHUFFLE(.,J). If $P=P_{R-1}P_{R-2}\cdots P_0$, according to the previous discussion regarding the action of UNSHUFFLE, we have that $A \bmod 2^R = P_{R-J-1}P_0\cdots P_{R-J+1}\cdots P_{R-1}$; from this observation, we readily obtain

$$A = L \cdot 2^R + (\tilde{P} \bmod 2^J) + (P \bmod 2^{R-J}) \cdot 2^J.$$

With regard to DLOOPOPER, we may reverse the sequence of operations performed in ALOOPOPER. Specifically, the data will be initially rearranged by applying OED; next, we have a sequence of operation-unshuffle pairs, and, at the end, the results are found in their correct order.
Specifically :

```
proc DLOOPOPER(L)
    OED(L);
    for J ← R-1 step -1 until J=0
        do foreach Q: 0≤Q<2<sup>R</sup>, BIT<sub>0</sub>(Q)=0
            pardo OPER(A,J; T[L·2<sup>R</sup>+Q], T[L·2<sup>R</sup>+Q+1])
            where A=L·2<sup>R</sup>+(Q mod 2<sup>J</sup>) + (Q mod 2<sup>R-J</sup>)/2<sup>J</sup>.
            odpar
        UNSHUFFLE(L,J)
    od
proc DLOOPOPER.
```

### 3.4) Sequential programs for each module

From the preceding global descriptions of `ASCEND` and `DESCEND`, it is straightforward to produce the sequential program of module(`ℓ`,p). The program `DMODULE(L,P)` for a given `DESCEND` algorithm is of the form : `HIGH SHEAVES(L,P)`; `LOWSHEAVES(L,P)`. Program `HIGH SHEAVES`, which processes sheaves k-1,...,r reads :

```
for T ← 1 step 1 until T=2<sup>R+1</sup>
do if 2<sup>R</sup> ≤ T+P < 2<sup>R+1</sup> then OPHIGH(L·2<sup>R</sup>+(P+T-1) mod 2<sup>R</sup>) fi;
    BSHIFT(L,P)
od
```

Program for `OPHIGH(M)` uses two registers `R1` and `R2` in each module and takes the form :

```
if BIT<sub>p</sub>(L)=0 then R2←Elink;
    (R1,R2)←function<sub>M,P+R</sub>(R1,R2);
    Elink←R2
else Elink←R1; R1←Elink
```

fi.

With these communication conventions, register `R1` always contains the module operand `T[ ]`.

Synchronization in the system must occur each time a communication link is used, so as to ensure that both the `get` and the `put` modules on each side of the link are ready for the transaction.
We omit the straightforward description of the programs for `SHIFT(L,P)`.
The program for `LOWSHEAVES` requires some care with respect to timing:

```pascal
proc LOWSHEAVES(L,P)
  for I ← R step -1 until I=2
  do UNSHUFFLE(I) od; comment end OED;
  for J ← R-1 step -1 until J=0
  do OPLOW(M,J) ; UNSHUFFLE(J+1)
    where M=L.2<sup>R</sup>+(P mod 2<sup>J</sup>) + (P̃ mod 2<sup>R-J</sup>).2<sup>J</sup>.
  od
  corp.
```

Here, `OPLOW(M,J)` performs:

```pascal
if BIT<0>(P)=0 then R2←Flink
  (R1,R2)←function<sub>M,J</sub>(R1,R2) ;
  Flink←R2
else Blink←R1 ; R1←Blink
```

As for `UNSHUFFLE(I)`:

```pascal
for T←1 step 1 until 2<sup>I-1</sup>
do if (T≤P<2<sup>I</sup>-T)∧(T=P mod 2)
  then R2←Flink ; Flink←R1 ; R1←R2
  elseif (T≤P<2<sup>I</sup>-T)∧(T≠P mod 2)
    then Blink←R1 ; R1←Blink
  fi
od.
```

The entire program for module (L,P) is of a very simple nature:
it basically counts up time T; at each step, it tests a simple logical condition involving T,L and P; depending on this tests, it performs an action
which can be of three types:

1) `get` operand ; compute function ; `put` operand

2) `put` operand ; `get` operand

3) do nothing.
The precise execution time of `ASCEND` (or `DESCEND`) on the CCC is given by the formula :

$$T = 4.2^R \cdot T_{CCC} + (R+2)^R T_{oper}$$
where $T_{CCC}$ is the time required for stepping up the control variable, testing it and performing one data exchange on some of the links ; $T_{oper}$ is the time required for computing $(U,V) \leftrightarrow \text{function}_{M,J}(U,V)$.
## 4. APPLICATIONS

In this section, we present applications of the CCC to the solution of interesting problems. Some applications -such as *bitonic merge* and *cyclic shift*-are directly within the ASCEND or DESCEND classes (simple algorithms) ; for these applications, all we have to do is specify OPER(M,J; U,V), since the individual programs for each module can be automatically derived from it, following the rules given in section 3. Other applications (such as *shuffle*, *unshuf*-le, *odd-even-distribution*, *odd-even-merge*, *fast-Fourier Transform*, *convolution*, *matrix transposition*) have programs which are cascades of a fixed (small) number of applications in the preceding class (*cascaded algorithms*); the corresponding algorithms on the CCC all run in $O(\log n)$ parallel steps.

We also have applications -such as *bitonic sort*, *odd-even sort*, and *calculation of symmetric functions*- for which the combining step of the two results of a recursive call is itself an algorithm in one of the two preceding categories. These algorithms, which we call *composite*, run in $O((\log n)^2)$ parallel steps on the CCC. We also mention that, in some instances, improvements can be found with respect to the general implementation scheme, especially in the LOOPOPER part of the program.

Finally, the CCC can also be used to solve a wide class of unrelated important problems, as we demonstrate with matrix multiplication.

### 4.1) Bitonic merge

The elegant algorithm for *bitonic merge* due to K. E. Batcher, [15] is ideally suited for implementation within the DESCEND class. Recall that a *bitonic sequence* $T[O:N]$ is either (1) increasing up to some point : $T[I] \le T[I+1]$ for $0 \le I <$, then decreasing from there on : $T[I] \ge T[I+1]$ for $P \le I < N$; or (2) a cyclic shift of a sequence as in (1).

If $T[O:N]$ is a bitonic sequence of even length (i.e., $N+1=2S$), the *bitonic split* operation.

foreach P: $0 \le P < S$ pardo COMPEXCHANGE(T[P], T[P+S]) odpar,
generates two bitonic sequences $T[O : S-1]$ and $T[S : 2S-1]$ such that each term of the former is no larger (resp. no smaller) than each term of the latter, provided COMPEXCHANGE(U,V) sets $U = \min(U,V)$ (resp. $max$).

Since the concatenation of two sequences $T[O : 2^{k-1}-1]$, $T[2^{k-1} : 2^k - 1]$ sorted in *opposite orders* forms a bitonic sequence, we can effectively merge these two sequences by performing successive bitonic splits on bitonic sequences
of lengths $2^k, 2^{k-1}, ..., 2$. This operation of bitonic merge can thus be implemented in the DESCEND class by specifying $OPER(A,J;U,V) = COMPEXCHANGE(U,V)$. In order to be able to use a bitonic merge for sorting, it will be convenient to define a type of bitonic merger which optionally realizes an increasing or a decreasing ordering. This new function, called ORIENTBITONICMERGE(L,I), is specified by the $OPER = ORIENTCOMPEXCHANGE$ where $ORIENTCOMPEXCHANGE(A,I; U,V)$ is

```pascal
if BIT_I(A)=0 then (U,V)=(min(U,V), max(U,V))
else (U,V)=(max(U,V), min(U,V))
fi
```

Of course, the running time of DLOOPOPER(L) can be improved with respect to the general format given in 3.3 by using a straightforward "transposition sort" discussed in section 4.4.

## 4.2) Radix-2 Fast Fourier Transform

The important FFT algorithm can be set in the ASCEND class. Let $\omega$ be a primitive root of unity of order $n=2^k$. The Fourier transform $<A_0, ..., A_{n-1}>$ of vector $<a_0, ..., a_{n-1}>$ is defined by

$$A_j = \sum_{0 \le i < 2^k} a_{2i} \omega^{ij} \quad \text{for } 0 \le j < 2^k.$$

A recursive description of the FFT algorithm readily follows from the algebraic identities :

for $0 \le j < 2^{k-1}$:
$$ 
\begin{cases}
A_j = \sum_{0 \le i < 2^{k-1}} a_{2i}(\omega^2)^{ij} + \omega^j \sum_{0 \le i < 2^{k-1}} a_{2i+1}(\omega^2)^{ij} ; \\
A_{j+2^{k-1}} = \sum_{0 \le i < 2^{k-1}} a_{2i}(\omega^2)^{ij} - \omega^j \sum_{0 \le i < 2^{k-1}} a_{2i+1}(\omega^2)^{ij}.
\end{cases}
$$

More concisely, we write

$A_j = U_j + \omega^j V_j$ and $A_{j+2^{k-1}} = U_j - \omega^j V_j$

where the U's and V's are respectively the Fourier Transforms, with respect to $\omega^2$, of the even subsequence $<a_0, a_2, ..., a_{2^{k-2}}>$ and the odd subsequence $<a_1, a_3, ..., a_{2^{k-1}}>$; we call $\omega^j$ the combining root powers.
The above relationships indicate that sequence $<a_0, ..., a_{n-1}>$ must be initially rearranged by means of the OED algorithm (given in 4.3). Once the desired reconfiguration has been achieved, we may proceed with the actual FFT computation, which is in the ASCEND class.

Its basic operation OPER(M,J;U,V) is specified by :

$$ (U,V) \leftarrow (U+\alpha V, U-\alpha V) \text{ where } \alpha = \omega^{M.2^{K-J}} $$

In order to complete the implementation of the FFT algorithm on the CCC, we must show that $\alpha$ can be computed efficiently at each step : a naive implementation computes $\alpha$ from $\omega$ by successive squaring, in time $O(k)$; since there are $k$ successive such computations within each module, the total time is $O(k^2)=O(\log^2 n)$. To reduce this time to $O(\log n)$, we consider the programs LOWSHEAVES and HIGHSHEAVES associated to module (L,P), by the transformations of section 3.

In HIGHSHEAVES, the successive values of the combining root power $\alpha$ are

$$ \omega(L.2^R+(P-I) \mod 2^R).2^{K-P-R} \quad \text{for } 0 \le I < 2^R $$

To keep the appropriate value of $\alpha$ in register $R\Omega$, we initialize $R\Omega \leftarrow \omega(L.2^R+P)$ and update $R\Omega \leftarrow R\Omega \times \omega^{2^{K-P-R}}$ every time $I$ is increased by one.

In LOWSHEAVES, the successive combining root powers are

$$ \alpha = \omega(L.2^R + (\tilde{P} \mod 2^I) + (P \mod 2^{R-I}) \cdot 2^I) 2^{K-I} \quad \text{for } 0 \le I < R $$

an expression obtained in section 3.3.

After simplification, we find

$$ \alpha = \omega^{\tilde{P}.2^{K-I}} = \omega_0^{2^{R-I}} \quad \text{with} \quad \omega_0 = \omega^{\tilde{P}.2^{K-R}} $$

Using R-I repeated squaring of $\omega_0$ allows to compute $\alpha$ in $O(r)$ steps; since there are $r$ successive operations of this type, timing for this part of the system is $O(R^2) = O(1\log n^2) = O(\log n)$.

### 4.3) Data rearrangements

Being able to efficiently permute the data on the CCC is obviously important for many applications. For example, the OED distribution is a necessary preliminary step to the FFT algorithm of the preceding section; note also that,
because of the duality indicated in section 3.1, a DESCEND version of FFT
could also be used; in this case, OED must be applied *after* rather than
*before* FFT.

### 4.3.1 Cyclic shifts

The procedure `SHIFT(T)` performs a cyclic forward shift of the array of
data `T[0:2^K-1]` by `T` positions. This procedure admits implementations in both
classes `ASCEND` and `DESCEND`.

The argument leading to a recursive description of `SHIFT` in the `ASCEND`
class is sketched in Figure 8 below.

<figure><img src="image_12.png" /></figure>

**Figure 8** A recursive method for shifting the network contents by `T` positions (case `0<T≤n/2`).

`SHIFT` is entirely specified by :

$$ \text{OPER}(M,J; U,V) = \begin{cases} 2^{J-1} \le (M \mod 2^{J-1}) + (T \mod 2^J) < 2^J & \text{if} \\ U \leftrightarrow V & \text{then} \\ \beta & \text{fi.} \end{cases} $$

Of course, program `LOOPOPER(L)` can be replaced by the more efficient version :

```
for I←1 step 1 until I=T mod 2^R
    do FSHIFT(L) od.
```

### 4.3.2 Shuffle, unshuffle and odd-even distribution

We have already considered `shuffle`, `unshuffle`, and `odd-even distribution` as rearrangements within a loop of the CCC. We shall now consider them at the level of the entire CCC.
We begin with "unshuffle", which is best explained by making reference to the $2^K \times K$ address array, defined by $T[M]=M$, that we have considered before (section 3.3). Define `REFLECT(S,I)` as the operation which takes the array $T[S:S+2^I-1]$ (where S is a multiple of $2^I$) and performs the following sequence of operations on it:

(i) the rows are permuted so that the leftmost (K-I) columns are left unaltered and the rightmost I columns give the standard gray-coding of the integers 0,1,...,2^I-1;

(ii) the rightmost I columns are cyclically shifted one position to the left among themselves (see figure 9 a,b,c).

<figure><img src="image_13.png" /><figcaption>Figure 9: Effect of `REFLECT(8,3)` on $T[8:15]$.</figcaption></figure>

Suppose now that the rows of $T[S:S+2^I-1]$ are initially permuted as specified by:

for each B : $S \le B < S + 2^I$ pardo if $BIT_0(B)=1$ then $T[B-1] \leftrightarrow T[B]$ fi od par
Obviously if we apply `REFLECT(S,I)` to this prearranged array, the resulting array is the mirror image with respect to its first row of that obtained by applying `REFLECT(S,I)` without any initial prearrangement (figure 9 c,e). Thus we have the following recursive definition of `REFLECT(S,I)` (see figure 10).

proc `REFLECT(S,I)`

1. for each B: $S+2^I-1 \le B < S+2^I$ par do if $BIT_0(B)=1$ then $T[B-1] \leftrightarrow T[B]$ fi od par

2. par do `REFLECT(S,I-1)`, `REFLECT(S+2^I-1, I-1)` od par

3. for each B: $S \le B < S+2^J$
par do if $BIT_{I-1}(B)=BIT_{I-2}(B)$ then $T[B] \leftrightarrow T[B-2^I-1]$ od par

corp `REFLECT`
<figure><img src="image_14.png" /><figcaption>**Figure 10:** Illustration of REFLECT(8,3) on T[8:15].</figcaption></figure>

In the above program, step 1 involves transpositions between adjacent modules; since it is followed by a recursive call, we can combine all steps 1 of these calls down to the deepest level of recursion, and obtain

OPER(M,0; U,V) = if $BIT_0(M)=1$) ∧ ($\sum_{1 \le J < K} BIT_J(M)=1$) then $T[M-1] \leftrightarrow T[M]$.

On the other hand, for $J \ne 0$ we have:

OPER(M,I; U,V) = if $BIT_I(M)=BIT_{I-1}(M)=1$ then $T[M] \leftrightarrow T[M-2^I]$.

Finally, we recall that the array representing UNSHUFFLE is given by the one-position column shift of the initial one. Thus, all we must do is transform the array obtained from REFLECT so as to leave the rightmost column unaltered and restore the binary coding from the Gray coding. It is easily verified that this is done by an ASCEND algorithm RESTORE specified by

OPER(M,I; U,V) = if $BIT_I(M)=BIT_{I-1}(M)=1$ then $T[M] \leftrightarrow T[M-2^{I-1}]$.

We conclude that UNSHUFFLE is of the form REFLECT; RESTORE and runs in time $O(\log n)$. An analogous treatment can be given for SHUFFLE, but will be omitted here.

We have already seen that OED can be achieved as a sequence of unshuffles (or of shuffles). However, rather than resorting to such a composite $O((\log n)^2)$ algorithm, we note that there are more direct $O(\log n)$ solutions, such as the following. For $0 \le B < 2^K$ define

LOW(B) = <BIT$_{\lfloor k/2 \rfloor - 1$}(B), ..., BIT$_0(B)$> and HIGH(B) = <BIT$_{k-1}$}(B), ..., BIT$_{\lfloor k/2 \rfloor}$}(B)>,
so that the binary representation of B is $BIN(B) = HIGH(B)LOW(B)$. We can regard
the sequence $T[0], ..., T[2^k-1]$ as made of consecutive blocks of size $2^{\lfloor k/2 \rfloor}$;
LOW(B) is then the address of B within the block, and HIGH(B) is the block address.
The address B where OED sends T[B] is at the HIGH(B) position within block LOW(B) (see figure 11).

By initially applying a right-cyclic shift to the contents of the block by

$$ \text{PSI}(B)=2^{\lfloor k/2 \rfloor} - \sum_{0 \le j < \lfloor k/2 \rfloor} 2^j \cdot \text{BIT}_{k-1-j} (B)=2^{\lfloor k/2 \rfloor} - (\tilde{B} \mod 2^{\lfloor k/2 \rfloor}) $$

positions, we correctly place the operands with respect to their final relative blocks addresses. We then perform a routing between blocks, processing the sheaves $\lfloor k/2 \rfloor, \lfloor k/2 \rfloor+1, \dots, k-1$. The final stage performs the inverse shifts of $2^{\lfloor k/2 \rfloor}$ - PSI(B) positions, within blocks.

T
0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
HIGH
0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1
0 0 0 0 1 1 1 1 0 0 0 0 1 1 1 1 1
LOW
0 0 1 1 0 0 1 1 0 0 1 1 0 0 1 1 1
0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1 0 1

T
0 1 2 3 6 7 4 5 11 8 9 10 13 14 15 12
0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1
0 0 0 0 1 1 1 1 0 0 0 0 1 1 1 1 1
0 0 1 1 1 1 0 0 1 0 0 1 0 1 1 0 0
0 1 0 1 0 1 0 1 1 0 1 0 1 0 1 0 0

T
0 8 4 12 6 14 2 10 13 1 9 5 11 7 15 3
0 1 0 1 0 1 0 1 1 0 1 0 1 0 1 0 1 0
0 0 1 1 1 1 0 0 1 0 0 1 0 1 1 0 0
0 0 0 0 1 1 1 1 0 0 0 0 1 1 1 1 1
0 0 0 0 0 0 0 1 1 1 1 1 1 1 1 1 1

T
0 8 4 12 2 10 6 14 1 9 5 13 3 11 7 15
shifts within blocks

<caption>Figure 11: The OED algorithm.</caption>

We can specify OED by the program : `SHIFT`; `ROUTING`; `UNSHUFFLE` where `SHIFT` corresponds to :

foreach P: $0 \le P < 2^{\lfloor k/2 \rfloor}$

pardo `SHIFT(P.2^{\lfloor k/2 \rfloor}, \lfloor k/2 \rfloor`, PSI(P.2^{\lfloor k/2 \rfloor})) odpar.
As for the ASCEND-program ROUTING, it is fully specified by

OPER(B,J) = if J >= ⌊k/2⌋ ∧ BIT<sub>K-J</sub>(B+PSI(S)) = 1
then EXCHANGE(B,J)
fi.

and its LOOPOPER program is void. Let us mention the interesting exercise of
pipelining the SHIFT computation with that of routing ; the idea to program the
initial SHIFTS in ASCEND and the final ones in DESCEND so that two successive
pipelining may smoothly take place.

### 4.3.3 Matrix transposition

Another rearrangement of inputs, useful in matrix algorithms, is that of
matrix transposition. Assume that T[0],...,T[2<sup>k</sup>-1] stores matrix A=[a<sub>i,j</sub>]<sup>k/2</sup> with
0≤i,j<2<sup>k/2</sup> (k even) in row major order, i.e. T[1+2<sup>k/2</sup>·j]=a<sub>i,j</sub> for 0≤i,J<2<sup>k/2</sup>.
In order to rearrange T in column major order, we perform an initial cyclic
shift of row (.,j) by j positions. A simple routing algorithm, using sheaves
k-1,k-2,...,k/2 then transfers a. to the appropriate column (figure 12).

Row major (0,0)(0,1)(0,2)(0,3) (1,0)(1,1)(1,2)(1,3) (2,0)(2,1)(2,2)(2,3) (3,0)(3,1)(3,2)(3,3)
Shifts (0,0)(0,1)(0,2)(0,3) (1,3)(1,0)(1,1)(1,2) (2,2)(2,3)(2,0)(2,1) (3,1)(3,2)(3,3)(3,0)
Routing (0,0)(1,0)(2,0)(3,0) (3,1)(0,1)(1,1)(2,1) (2,2)(3,2)(0,2)(1,2) (1,3)(2,3)(3,3)(0,3)
Shifts (0,0)(1,0)(2,0)(3,0) (0,1)(1,1)(2,1)(3,1) (0,2)(1,2)(2,2)(3,2) (0,3)(1,3)(2,3)(3,3)

<caption>Figure 12 : In situ matrix transposition.</caption>

proc TRANSPOSE
  foreach J: 0≤J<2^K/2 pardo SHIFT(2^K/2·J,K/2,J) odpar;
  ROUTE(0,K-1);
  foreach J: 0≤J<2^K/2 pardo SHIFT(2^K/2·J,K/2,-J) odpar
  corp TRANSPOSE.

Here, the DESCEND-algorithm ROUTE is specified by OPER(M,J; U,V):

if (J ≥ K/2) ∧ BIT<sub>K-J-1</sub>(M+M div 2<sup>K/2</sup>) = 1 then U ↔ V fi.

This algorithm is very close to the second OED algorithm presented earlier :
indeed, if the binary representation of M is HIGH-LOW, OED exchanges M with
LOW-HIGH while TRANSPOSE exchanges it with LOW-HIGH.
### 4.4) Odd-even-merge

Odd-even merge is a procedure, also introduced by Batcher [15], to merge in parallel two sorted sequences. It derives its name from the feature that odd-indexed and even-indexed keys in either sequence are treated separately. The traditional description is as follows:

(1) the two sequences are concatenated and unshuffled (which achieves the odd-even split)

(2) the two resulting pairs of sequences are merged recursively, and then

(3) combined into one sorted sequence by a single parallel step.

Specifically, if A[1:T] and B[1:T] are the two sorted sequences obtained in (2), step (3) consists in

```pascal
foreach J:2≤J≤T   pardo COMPEXCHANGE(A[J],B[J-1]) odpar
```

It is easily realized that the succession of odd-even splits (i.e., of unshuffles) corresponds exactly to the GOED procedure described in section 3. Thus, we assume that the two initial sorted sequences be subjected to a preparatory GOED.

Subsequently, after performing a single comparison-exchange of each even-indexed key with the following odd-indexed one, we must combine, according to (3) above, adjacent strings of two keys, then of four, and so on in a standard hierarchical fashion. The ensuing algorithm, OECOMBINE (O,K) is obviously in the ASCEND class and will be entirely defined by its OPER specification.

However, in order to improve efficiency, we may want to replace the standard ALOOPOPER by an "ad hoc" sorting procedure which is naturally suited to linearly connected arrays as the CCC loops. Additionally in order to correctly implement the 'sheaf' portion of OECOMBINE, it is essential to realize the correct positioning of the initial orderings in each loop (according to (3) above). Specifically, we say that loop *e* has *offset* *v*(*e*) if the first element of the sorted sequence stored in *e* occurs in module (*e*,v(*e*)). It is easy to realize that v(*e*)=v(*e*)+1 if and only if the binary expansions of *e* and *e'* differ in exactly one position.

Thus, choosing v(0)=0, we obtain

$$v(e) = \sum_{0 \le j < k - r} \text{BIT}_j(e).$$
The initial ordering in loop *e* can therefore be achieved by the following
"ad hoc" loop procedure (which is not classifiable in either class) :

```
proc OFFSETTRANSPSORT(L)
    for I<1 step 1 until I=2^R-1
        do foreach P:P<>v(L), P=2Q+I(mod 2), 0≤Q<2^R-1
            pardo COMPEXCHANGE(T[P] T[(P-1)mod 2]^R, odpar
        od
    corp OFFSETTRANSPSORT.
```

Notice that by "breaking" the loop between module $e.2^I+v(e)-1$ and $e.2^I+v(e)$ we realize at the same time a sorted arrangement and a cyclic forward shift of $v(e)$ positions.

The OPER specification for OECOMBINE (O,K) is OPER(M,J;U,V)=COMPEXCHANGE(U,V).

In conclusion we have, OEMERGE (O,K)=GOED(O,K) ; OECOMBINE (O,K).

## 4.5 Sorting algorithms

The previously described merging routines may obviously be used to implement
sorting algorithms.

At first we consider a sorting procedure based on bitonic merger. We call
ORIENTLOOPBITONIC(L) the DLOOPOPER corresponding to ORIENTBITONICMERGE(O,K).
We immediately obtain :

```
proc SORT(P,I)
    if I>0 then
        pardo SORT(P,I-1), SORT(P+2^(I-1),I,1) odpar ;
        ORIENTBITONICMERGE(P,I)
    fi
corp SORT
```

An analogous algorithm is based on OEMERGE. Both algorithms run in time $O((\log n)^2)$.

## 4.6 Convolution and calculation of symmetric functions

Given two sequences $U[S:S+2^{I-1}-1]$ and $U[S+2^I:S+2^{I-1}]$, CONVOLVE(S,I) computes
their convolution and has the following format
**proc CONVOLVE(S,I)**

$T[S:S+2^{I-1}-1] \leftarrow U[S:S+2^{I-1}-1], T[S+2^{I-1}:S+2^{I-1}] \leftarrow 0$

$M[S:S+2^{I-1}] \leftarrow \text{FFT}(S,I)$

$T[S:S+2^{I-1}] \leftarrow U[S+2^{I-1}:S+2^{I-1}], T[S+2^{I-1}:S+2^{I-1}] \leftarrow 0$

$T[S:S+2^{I-1}] \leftarrow \text{FFT}(S,I)$

$\underline{\text{foreach } B: S \le B < S+2^{I} \text{ pardo } T[B] \leftarrow T[B] \times M[B] \text{ odpar}}$

**corp CONVOLVE**

Obviously, CONVOLVE is a cascaded algorithm and runs in time $O(\log n)$. As to the calculation of symmetric functions we have

**proc SYMM(P,I)**

if $I > 0$ then

$\underline{\text{pardo SYMM}(P, I-1), \text{SYMM}(P+2^{I-1}, I-1) \text{ odpar}}$

**corp SYMM**

which has exactly the same structure as SORT above.

## 4.7 Matrix multiplication

To compute the matrix product $C=A \times B$ of two $n \times n$ matrices ($\omega \log n = 2^i$), we must obviously first store

$$A = \begin{pmatrix} A_0 & & & \\ & \ddots & & \\ & & A_{n-1} & \\ & & & A_n \end{pmatrix}$$ in row major order, and $B = (B_0 \cdots B_{n-1})$ in column major order

(see section 4.3.3). Assuming the CCC large enough, i.e. $2^k \ge n^3$, we copy $A$ and $B$ into the pattern :

$$A_0 B_0 A_{n-1} B_0 \cdots A_0 B_j \cdots A_{n-1} B_{n-1}$$

All this can be achieved with simple-minded algorithms, in time $O(\log n)$ on the CCC. Each of the scalar products.

$c_{i,j} = A_i \cdot B_j = \sum_k a_{i,k} b_{k,j}$ is computed in parallel, within $O(\log n)$ additional time units. The results $c_{i,j}$ are then regrouped, according to the output format (say, row major), performing one last rearrangement of the CCC very similar to TRANSPOSE described above. Although the details of this algorithm are a bit tedious to describe, it should be clear that *matrix multiplication* can be *computed on the CCC in time O(\log n)*.
Of course, if the CCC is not large enough, say $2^k=2 \cdot n^2$, we do not have
space for performing the initial copying of A and B. The best we can hope for
in that case is an $O(n)$ time bound since (naive) matrix multiplication requires
$O(n^3)$ arithmetic operations, only $O(n^2)$ of which can simultaneously take place
at a given time. The linear time bound $O(n)$ can be achieved on the CCC, in a
rather straightforward manner, but we shall not attempt to describe it here.
## 5 - LIMITED PARALLELISM

So far, we have assumed that the size n of the CCC was tailored to the application. In many cases, it is more realistic to assume that once the system is constructed, one may wish to use it for solving problems where the number N of inputs is not necessarily equal to n, but typically N may be larger. To cope with this kind of situation, we suggest to let each module of the CCC be a full fledged microprocessor endowed with a private RAM memory.

Assuming for simplicity that $N=2^n$, with $t=2^q$ integer, we require that the RAM memory of each module be of size t.

Let us denote by $T[M,0]$, $T[M,1]$, ..., $T[M,2^q-1]$ the private memory locations of module M. For specificity, we describe the treatment of a DESCEND program; the other applications are straightforward modifications.

The input $a_0, ..., a_N$ is divided into consecutive blocks of size $t=2^q$, each block being stored within a module of the CCC; specifically, $T[M,J]$ contains $a_{2^q \cdot m + j}$ for $0 \le j < 2^q$.

Program for module M(L,P) (see section 3.4) has the form HIGHSHEAVES(L,P); LOWSHEAVES(L,P) ; TERMINAL(L,P).

In the programs for HIGHSHEAVES and LOWSHEAVES given in section 3.4, we merely change OPHIGH and OPLOW respectively to :

```
proc NEWOPHIGH(M)
    for K<0 step 1 until K=2^q-1
    do R1←T[M+K]; OPHIGH(M+K); T[M+K] ← R1 od
```
```
corp NEWOPHIGH ;
```

Similarly, NEWOPLOW is obtained by replacing HIGH by LOW throughout.

The program TERMINAL(L,P) which terminates the algorithm within module M has the obvious form :

```
proc TERMINAL(L,P)
    M 2^Q. (P+2^R.L) ;
    for J←Q-1 step -1 until J=0
    do for I←0 step, until I=2^q-1
        do if BIT_J(I)=0 then OPER(M+I,J; T[M,I], T[M,I+2^J]) fi od
        od
    od
    corp TERMINAL.
```
It should be clear by now that all of the algorithms described in section 4 can be applied here, and we have a very efficient and economical means for connecting microprocessors. A direct analysis shows that, on a CCC comprising n processors, each processor having memory $\frac{N}{n}$, we can process N inputs in time $O(\frac{N}{n} \log n)$ for algorithms in the classes ASCEND or DESCEND.
6 - CONCLUSION

In this paper, we have proposed a structure which can be used for direct hardware implementation of specific useful algorithms, or, as suggested in section 5, as a general purpose parallel processing system.

We have limited our presentation to the algorithms which cannot be processed efficiently under simpler connections ; it should be clear by now that, in our general purpose system, algorithms for matrices and graphs problems can also be coded as efficiently as on any other known system. The proposed architecture is thus both quite powerful and within the possibilities of present technology. It raises both practical problems, regarding the implementation or physical realization, and theoretical problems : can we obtain complexity results (e.g. can sorting be performed in time O(logn) on our network ?) in a machine model that would be general enough to include our scheme ? Can we specify a high level programming language for parallel computation, which would be natural and simple to use, and nevertheless automatically compilable on devices such as ours ? Can all of the known efficient parallel algorithms be implemented on the CCC without significant loss of efficiency ?
7 - REFERENCES

[1] A.V. Aho, J.E. Hopcroft and J.D. Ullman,
*The Analysis and Design of Computer Algorithms*
Addison-Wesley, Reading, Mass., 1974.

[2] D. Heller,
"A Survey of parallel algorithms in numerical linear algebra"
Dept. of Comp. Sci. Carnegie-Mellon University,
Pittsburgh, Pa., Feb. 1976.

[3] D.S. Hirschberg,
"Fast parallel sorting algorithms"
Tech. Rep, Dept. of Electr. Eng.,
Rice University, Houston, Tx., Jan. 1977.

[4] L.G. Valiant,
"Parallelism in Comparison Problems"
SIAM Journal on Computing, Vol 4, 3, pp. 348-355, Sept. 1975.

[5] F.P. Preparata,
"New parallel sorting schemes"
IEEE Transaction on Computers,
Vol. C-27, No. 7, pp. 669-673, July 1978.

[6] W.M. Gentleman,
"Some complexity result for matrix computation on parallel processors"
Jour. of the ACM, 25,1,pp. 112-115, January 1978.

[7] G.H. Barnes et al,
"The ILLIAC IV computer"
IEEE Transactions on Computers, C-17, pp. 746-757, 1968.
[8] C.D. Thompson and H.T. Kung,
"Sorting on a mesh connected computer"
Proc. ACM-SIGACT Symp. on Theory of Computing,
Hersley, Pa., pp. 58-64, May 1976.

[9] H.T. Kung and C.E. Leiserson,
"Algorithms for VLSI processor arrays"
Symposium on Sparse Matrix Computations, Knoxville, Tenn., Nov. 1978.

[10] D. Nassimi and S. Sahni,
"Bitonic Sort on a Mesh-Connected Parallel Computer"
IEEE Transaction on Computers, Vol C-28, No. 1, pp. 2-7, Jan. 1979.

[11] L.J. Guibas, H.T. Kung, and C.D. Thompson,
"Direct VLSI implementation of Combinatorial algorithms"
Research Report, Department of Computer Science,
Carnegie-Mellon University, Pittsburgh, Penn., March 1979.

[12] K.N. Levitt and W.H. Kautz,
"Cellular arrays for the solution of graph problems"
Comm. of the ACM, vol 15, No. 9, pp. 789-801, 1972.

[13] M.C. Pease,
"The indirect binary n-cube microprocessor array"
IEEE Trans. on Computers, Vol. C-26, No. 5, pp. 458-473, May 1977.

[14] A.M. Mead and L.A. Convay,
Introduction to VLSI Systems
Textbook in preparation (1978)

[15] K.E. Batcher,
"Sorting networks and their applications"
Proc. AFIPS Spring Joint Computer Conference, Vol 32, pp. 307-314, April 1968.

[16] D.E. Knuth,
The Art of Computer Programming
Volume 3, Sorting and Searching, Addison-Wesley,
Reading, Mass., 1973.
[17] H.S. Stone,
"Parallel processing with the perfect shuffle"
*IEEE Transactions on Computers*
Vol. C-20, pp. 153-161, 1971.

[18] C.D. Thompson,
"Area-Time Complexity for VLSI",
Proc of the 11th Annual ACM Symp. of Theory of Computing,
pp. 81-88, 1979.

[19] C.D. Thompson,
"A complexity theory for VLSI",
Ph.D. Thesis, Carnegie Mellon University, Dept of Computer Science, 1980.

[20] M.C. Pease,
"An adaptation of the fast Fourier transform for parallel processing",
*Journal of the ACM*, 15(2), pp 252-264, April 1968.