\documentstyle[12pt,fullpage]{article}
%\renewcommand{\baselinestretch}{1.4}
\title{Uniform Approach for Solving \\Some Classical Problems on a Linear Array}

\author{David R. O'Hallaron\thanks{An earlier version of this paper was presented at the International Conference on Parallel Processing, St. Charles, IL, August, 1989. This work was supported in part by the GE Aerospace Business Group, while the author was with the GE Research and Development Center, Schenectady, NY, 12301. The author's current address is School of Computer Science, Carnegie Mellon University, Pittsburgh, PA, 15213.}}

\date{}

\begin{document}
\baselineskip=18pt
\bibliographystyle{ieeetr}
\maketitle
\begin{abstract} 
We show that the same basic recurrence can be used to solve a number
of classical problems from linear algebra and graph theory, including
instances of the algebraic path problem, matrix multiplication, matrix
triangularization, and matrix transpose, and we give a simple mapping
of the recurrence onto a unidirectional linear array. Qualitative
advantages to programming linear arrays using this approach include
uniformity of design, simplicity of programming, and scalability to
larger problems.  The major disadvantage is that the resulting
algorithms are not necessarily optimal.
\end{abstract}


\newpage
\section*{Keywords} 
parallel algorithm, linear array, matrix computation, linear algebra,
graph theory.
\newpage

\section{Introduction}
In this paper, we describe a technique for solving a number of
classical problems on a unidirectional linear array using essentially
the same parallel algorithm.

First, we demonstrate that some interesting problems from linear
algebra and graph theory, including instances of the algebraic path
problem, matrix multiplication, matrix triangularization, and matrix
transpose, can be solved using the same recurrence.  Informally, at
the $kth$ step in the recurrence, transformation $T^{(k)}$ is applied
to square matrix $A^{(k-1)}$.  The resulting matrix is then permuted
by circularly shifting its columns one position towards the left and
its rows one position towards the top.  The permuted matrix $A^{(k)}$
becomes the input to the $(k+1)st$ step.

Next, we give a simple mapping of the recurrence onto a linear array,
where the $kth$ cell inputs matrix $A^{(k-1)}$ from the left, applies
transformation $T^{(k)}$, permutes the resulting matrix, and sends the
permuted matrix $A^{(k)}$ to the right, all on an element-by-element
basis.

This approach (when applicable) offers several qualitative advantages.
First, it can simplify the design process. If we are able to cast a
problem into the desired recurrence, the mapping onto the linear array
is straightforward and mechanical. Second, it can simplify the
programming process because each cell executes the same code. Third, it
is easily generalized to handle large problems on small arrays.
The major disadvantage of the approach is that the resulting
algorithms are not necessarily optimal.

The work builds on previous work by S.Y. Kung \cite{sykung87} and P.S.
Tseng \cite{pstseng87}, each of whom used permutations to design
parallel algorithms for instances of the algebraic path problem, Kung
using a 2-D array, Tseng using a linear array. Our contribution is to
show that a similar technique is useful on a larger class of problems.
These problems all share the same localized and regular dependency
structure that is exploited by systolic algorithms \cite{kung82a},
wavefront algorithms \cite{sykung81}, and automatic synthesis
algorithms \cite{quinton84,rao85,moldovan86,leekedem88,leekedem90}.

Section 2 describes the recurrence and gives numerous examples of its
utility. Section 3 describes a simple mapping of the recurrence onto a
linear array. The paper concludes with a discussion of the results.

\section{Recurrence}
A surprising number of classical problems from graph theory and linear
algebra can be solved using the same simple recurrence. We define the
recurrence and then give some examples of its usefulness.

\subsection{Definitions}
Given an $n \times n$ matrix $A = \{a_{ij}\}$, we want to compute the
$n \times n$ matrix, $B = \{b_{ij}\}$, using a recurrence of the form:
\begin{eqnarray}
A^{(0)} & = & A \nonumber\\
A^{(k)} & = & C T^{(k)} (A^{(k-1)}) C^{-1} \label{eqn:problem},
~~k = 1, \ldots ,m\\
B &=& A^{(m)} \nonumber
\end{eqnarray}
The symbol $C$ denotes an $n \times n$ permutation matrix, 
\begin{displaymath}
C = 
\left[
\begin{array}{ccccc}
e_n & e_1 & e_2 & \cdots & e_{n-1}
\end{array}
\right]
\end{displaymath}
where $e_i$ is an $n \times 1$ column vector such that the only
non-zero element is a $1$ in the $ith$ position.  Notice that
\begin{displaymath}
C^{-1} = C^{T} = 
\left[
\begin{array}{ccccc}
e_2 & e_3 & \cdots & e_{n} & e_1
\end{array}
\right].
\end{displaymath}
Multiplying a matrix $A$ on the left by $C$ and on the right by its
inverse $C^{-1}$ has the effect of circularly shifting each column of
$A$ one position towards the left and each row of $A$ one position
towards the top. For example,
\begin{displaymath}
\left[
\begin{array}{rrr}
5&6&4\\
8&9&7\\
2&3&1
\end{array}
\right]
= C
\left[
\begin{array}{rrr}
1&2&3\\
4&5&6\\
7&8&9
\end{array}
\right]
C^{-1}.
\end{displaymath}
S.Y. Kung \cite{sykung87} and P. S. Tseng \cite{pstseng87} have used a
similar permutation to design parallel algorithms for instances of the
algebraic path problem, Kung in the context of an optimal algorithm
for 2-D arrays, and Tseng in the context of an implementation of the
shortest paths problem on the Warp computer \cite{annaratone87e}.  As
we shall see, this type of permutation is useful for solving other
kinds of problems as well.

The symbol $T^{(k)}$ denotes a transformation applied to $A^{(k-1)}$
such that
\begin{displaymath}
\hat{A}^{(k)} = \{\hat{a}_{ij}^{(k)}\} = T^{(k)}(A^{(k-1)})
\end{displaymath}
where $\hat{A}^{(k)}$ is an $n \times n$ matrix. Let
\begin{equation}
S_{ij} = 
\left\{
(p,q) ~|~ \mbox{$p \leq \max(i, i^{\prime})$ and $q \leq \max(j, j^{\prime})$}
\right\}
\label{eqn:set1}
\end{equation}
be a set of {\em index pairs}, where for some index $r$,
\begin{displaymath}
r^{\prime} = 
\left\{
\begin{array}{ll}
n&\mbox{if $r = 1$}\\
r-1&\mbox{otherwise.}
\end{array}
\right.
\end{displaymath}
Let $(p,q)^{(k-1)}$ denote the index pair associated with
$a_{pq}^{(k-1)}.$ Then the only restriction on transformation
$T^{(k)}$ is that each $\hat{a}_{ij}^{(k)}$ must be computable using
only elements $a_{pq}^{(k-1)}$ such that $(p,q)^{(k-1)} \in S_{ij}$.

Set $S_{ij}$ defines a restriction on transformation $T^{(k)}$ that
limits the collection of elements of $A^{(k-1)}$ that can be used to
compute the element $\hat{a}_{ij}^{(k)}$. As we will see in the
following sections, the restriction on transformation $T^{(k)}$
defined by $S_{ij}$ is weak enough to allow a significant number of
problems from linear algebra and graph theory to be solved using a
recurrence in the form of Equation (\ref{eqn:problem}), yet strong
enough to admit a uniform parallel solution for each of these problems
on a linear array. 

The set $S_{ij}$ can be interpreted visually using the graph in Figure
\ref{fig:pgraph}, for the case where $n=4$.  Unfilled vertices in the
graph correspond to index pairs of $A^{(k-1)}$. Filled vertices in the
graph correspond to index pairs of $A^{(k)}$. If we let $a
\stackrel{\ast}{\rightarrow} b$ denote the existence of a path from
vertex $a$ to vertex $b$, then we can interpret set $S_{ij}$ as
\begin{equation}
S_{ij} = 
\left\{
(p,q)^{(k-1)} ~|~ (p,q)^{(k-1)} \stackrel{\ast}{\rightarrow}
(i^{\prime} , j^{\prime})^{(k)}
\right\}.
\label{eqn:set2}
\end{equation}
For example, using either Equation (\ref{eqn:set1}) or
(\ref{eqn:set2}), we find that for $n=4$,
\begin{eqnarray*}
S_{11} &=& \{(1,1), \ldots (1,4), (2,1), \ldots, (4,3), (4,4)\}\\
S_{12} &=& \{(1,1), \ldots (4,1), (1,2), \ldots, (4,2)\}\\
S_{22} &=& \{(1,1), (1,2), (2,1), (2,2)\}
\end{eqnarray*}
\subsection{Examples}
This section gives examples of some problems --- in particular,
instances of the algebraic path problem, matrix multiplication, matrix
triangularization, and matrix transpose --- that can be solved using the
recurrence in Equation (\ref{eqn:problem}). 

\subsubsection{Instances of the Algebraic Path Problem}
A number of problems, including the {\em shortest paths}, {\em
transitive closure}, and {\em minimum spanning tree} problems from
graph theory, and the {\em matrix inverse} problem from linear
algebra, have been shown to be instances of the general problem known
as the {\em algebraic path problem} (APP). Formal definitions of the
APP can be found in \cite{sykung87,rote85}.  All instances of the APP
can be solved using the recurrence in Equation (\ref{eqn:problem}).
We give three examples: shortest paths, minimum spanning tree, and
matrix inverse.

\paragraph{\bf Shortest Paths:}
We are given a directed graph $G$ with $n$ vertices and labeled with
non-negative edge weights.  The length of a path $p$ is defined to be
the sum of the edge weights along $p$.  The problem is to determine
the length of the shortest path between all pairs of vertices $a$ and
$b$.  In Equation (\ref{eqn:problem}), let $A$ be the adjacency matrix
of a graph $G$ containing $m=n$ vertices and let
\begin{displaymath}
\hat{a}_{ij}^{(k)} =
\min(a_{ij}^{(k-1)}, a_{i1}^{(k-1)} + a_{1j}^{(k-1)}), ~~1 \leq i,j \leq n.
\end{displaymath}
Then $a_{ij}^{(n)}$ is the shortest path from vertex
$v_i$ to vertex $v_j$ in $G$. This formulation is equivalent to
solutions presented in \cite{sykung87,rote85,floyd62}. 
For example, consider the graph in Figure \ref{fig:exgraph} with
adjacency matrix 
\begin{displaymath}
A = A^{(0)} =
\left[
\begin{array}{rrrr}
0&2&1&\infty\\
2&0&\infty&3\\
1&\infty&0&7\\
\infty&3&7&0
\end{array}
\right]
\end{displaymath}
Applying Equation (\ref{eqn:problem}) gives
\begin{displaymath}
A^{(1)} =
\left[
\begin{array}{rrrr}
0&3&3&2\\
3&0&7&1\\
3&7&0&\infty\\
2&1&\infty&0
\end{array}
\right]
\end{displaymath}
\begin{displaymath}
A^{(2)} = 
\left[
\begin{array}{rrrr}
0&6&1&3\\
6&0&5&3\\
1&5&0&2\\
3&3&2&0
\end{array}
\right]
\end{displaymath}
\begin{displaymath}
A^{(3)} =
\left[
\begin{array}{rrrr}
0&5&3&6\\
5&0&2&1\\
3&2&0&3\\
6&1&3&0
\end{array}
\right]
\end{displaymath}
\begin{displaymath}
A^{(4)} =
\left[
\begin{array}{rrrr}
0&2&1&5\\
2&0&3&3\\
1&3&0&6\\
5&3&6&0
\end{array}
\right]
\end{displaymath}
The reader can verify by inspection that $A^{(4)}$ contains the
shortest paths between all pairs of vertices.

\paragraph{\bf Minimum Spanning Tree:}
We are given an undirected graph $G$ with $n$ vertices, labeled with
distinct non-negative edge weights. A spanning tree of $G$ is a subgraph
of $G$ that contains $n$ vertices and $n-1$ edges. The size of a spanning
tree is the sum of its edge weights. The problem is to find the smallest
spanning tree of $G$.  
In Equation (\ref{eqn:problem}) let $A^{(0)} = A$ be the adjacency
matrix of a graph $G$ containing $m=n$ vertices and let
\begin{displaymath}
\hat{a}_{ij}^{(k)} =
\min(a_{ij}^{(k-1)}, ~\max(a_{i1}^{(k-1)}, a_{1j}^{(k-1)})), ~~
1 \leq i,j \leq n.
\end{displaymath}
Then by a result in \cite{maggs88}, edge $(i,j)$ is in the
minimum spanning tree of $G$ if and only if $a_{ij}^{(0)} =
a_{ij}^{(n)}$. 
For example, consider the graph in Figure
\ref{fig:exgraph}. Applying Equation (\ref{eqn:problem}) gives
\begin{displaymath}
A^{(4)} =
\left[
\begin{array}{rrrr}
0&2&1&3\\
2&0&2&3\\
1&2&0&3\\
3&3&3&0
\end{array}
\right]
\end{displaymath}
Comparing the elements of $A^{(4)}$ with the elements of $A^{(0)}$,
the reader will find that the minimum spanning tree of the graph in
Figure \ref{fig:exgraph} consists of the edges
$(v_{1},v_{2})$, $(v_{1},v_{3})$, and $(v_{2},v_{4})$.

\paragraph{\bf Matrix Inverse:}
We are given a nonsingular $n \times n$ matrix $X$ with non--zero
diagonal elements. The problem is to find (without pivoting) the
unique matrix $X^{-1}$ such that $XX^{-1} = X^{-1}X = I$. In Equation
(\ref{eqn:problem}), let $A = I - X$. Let $m=n$ and let
\begin{displaymath}
\hat{a}_{ij}^{(k)} =
\left\{
\begin{array}{ll}
{(1 - a_{11}^{(k-1)})}^{-1}&\mbox{if $i=1$ and $j=1$}\\*[.1in]
a_{1j}^{(k-1)}{\hat{a}_{11}^{(k)}}&\mbox{if $i=1$ and $j>1$}\\*[.1in]
a_{i1}^{(k-1)}{\hat{a}_{11}^{(k)}}&\mbox{if $i>1$ and $j=1$}\\*[.1in]
a_{ij}^{(k-1)} + \hat{a}_{i1}^{(k)}a_{1j}^{(k-1)}&\mbox{otherwise.}
\end{array}
\right.
\end{displaymath}
Then $X^{-1} = A^{(n)}$. Note that $\hat{a}_{11}^{(k)}$ is always
defined because of the assumption of non--zero diagonal elements in
$X$. This is equivalent to the formulation given in
\cite{robert86}. For example, consider the problem of finding the
inverse of
\begin{displaymath}
X = 
\left[
\begin{array}{rrr}
1&1&1\\
1&2&3\\
0&1&1
\end{array}
\right]
\end{displaymath}
We let 
\begin{displaymath}
A = A^{(0)} = I-X =
\left[
\begin{array}{rrr}
0&-1&-1\\
-1&-1&-3\\
0&-1&0
\end{array}
\right]
\end{displaymath}
Applying Equation (\ref{eqn:problem}) gives
\begin{displaymath}
A^{(3)} = 
\left[
\begin{array}{rrr}
1&0&-1\\
1&-1&2\\
-1&1&-1
\end{array}
\right]
\end{displaymath}
The reader can verify that $A^{(3)} = X^{-1}$ by using the identity
$XX^{-1} = I$.

\subsubsection{Matrix Multiplication}
We are given $m \times m$ matrices $X$ and $Y$. The problem
is to compute the $m \times m$ matrix $Z = XY$.
Using a technique similar to one described in \cite{aho74}, let $A$ in
Equation (\ref{eqn:problem}) be the $n \times n$ matrix
\begin{displaymath}
A = A^{(0)} = 
\left[
\begin{array}{rr}
0&Y\\
X&0
\end{array}
\right]
\end{displaymath}
where $n = 2m$. Let
\begin{displaymath}
\hat{a}_{ij}^{(k)} = 
a_{ij}^{(k-1)} + a_{i1}^{(k-1)} a_{1j}^{(k-1)}, ~~1 \leq i,j \leq n
\end{displaymath}
Then 
\begin{displaymath}
A^{(m)} =
\left[
\begin{array}{rr}
Z&X\\
Y&0
\end{array}
\right].
\end{displaymath}
For example, suppose that 
\begin{displaymath}
X = 
\left[
\begin{array}{rr}
1&2\\
-1&3
\end{array}
\right]
\end{displaymath}
and
\begin{displaymath}
Y = 
\left[
\begin{array}{rr}
2&1\\
0&1
\end{array}
\right]
\end{displaymath}
Then
\begin{displaymath}
A = A^{(0)} =
\left[
\begin{array}{rrrr}
0&0&2&1\\
0&0&0&1\\
1&2&0&0\\
-1&3&0&0
\end{array}
\right]
\end{displaymath}
and applying Equation (\ref{eqn:problem}) gives
\begin{displaymath}
A^{(1)} =
\left[
\begin{array}{rrrr}
0&0&1&0\\
2&2&1&1\\
3&-2&-1&-1\\
0&2&1&0
\end{array}
\right]
\end{displaymath}
\begin{displaymath}
A^{(2)} =
\left[
\begin{array}{rrrr}
2&3&1&2\\
-2&2&-1&3\\
2&1&0&0\\
0&1&0&0
\end{array}
\right]
\end{displaymath}
The reader can verify by inspection that
\begin{displaymath}
A^{(2)} = 
\left[
\begin{array}{rr}
Z&X\\
Y&0
\end{array}
\right]
\end{displaymath}

\subsubsection{Matrix Triangularization}
We are given an $n \times n$ matrix $A$. The problem is
to compute the $n \times n$ matrix product
\begin{equation}
B = A \Theta
\label{eqn:givens1}
\end{equation}
where $\Theta$ is any $n \times n$ orthogonal matrix that
triangularizes $A$. This important computation arises frequently in
estimation applications such as Kalman filtering, where a factored
form of a covariance matrix is updated in-place
\cite{jover86,baheti90}.  A popular technique for
computing $B$ is to use orthogonal rotations to zero the elements of
the upper triangle of matrix $A$\cite{golub83}.  An orthogonal
rotation is an $n \times n$ matrix of the form
\begin{displaymath}
\Theta_{ij} = 
\left[
\begin{array}{ccccc}
1&\vdots&&\vdots&\\
\cdots&r_{ii}&\cdots&r_{ij}&\cdots\\
&\vdots&1&\vdots&\\
\cdots&r_{ji}&\cdots&r_{jj}&\cdots\\
&\vdots&&\vdots&1\\
\end{array}
\right]
\end{displaymath}
where 
\begin{eqnarray*}
r_{ii} &=& \cos(\theta)\\
r_{ij} &=& \sin(\theta)\\
r_{ji} &=& \sin(\theta)\\
r_{jj} &=& -\cos(\theta)
\end{eqnarray*}
and where 
\begin{eqnarray*}
\theta &=& \arctan(a_{ii}/a_{ij})\\
\cos(\theta) &=& a_{ii} / (a_{ii}^2 + a_{ij}^2 )^{1/2}\\
\sin(\theta) &=& a_{ij} / (a_{ii}^2 + a_{ij}^2 )^{1/2}.
\end{eqnarray*}
Postmultiplying $A$ by $\Theta_{ij}$ zeros element $a_{ij}$. Since the
product of two orthogonal matrices is orthogonal, Equation 
(\ref{eqn:givens1}) can be computed by postmultiplying $A$ by 
the appropriate orthogonal rotations:
\begin{equation}
B = (( \cdots ((A \Theta_{12})\Theta_{13}) \cdots ) \Theta_{n-1,n})
\label{eqn:givens2}
\end{equation}
Equation (\ref{eqn:givens2}) can be computed using the recurrence in
Equation (\ref{eqn:problem}). Let $A$ be the $n \times n$ matrix. Let
$m = n$ and let
\begin{displaymath}
<p,q>_r ~\equiv~
\left[ 
\begin{array}{ccc}
a_{p1}^{(k-1)} & \cdots & a_{pr}^{(k-1)} 
\end{array}
\right]
\left[
\begin{array}{c}
a_{q1}^{(k-1)}\\
\vdots\\
a_{qr}^{(k-1)}
\end{array}
\right]
\end{displaymath}
In words, $<p,q>_r$ is an inner product formed by multiplying the 
row vector formed from the first $r$ elements of row $p$ by the 
transpose of the row vector formed from the first $r$ elements of
row $q$. Also, let
\begin{displaymath}
n(i,j,k) ~\equiv~ {<1,i>_{j-1}a_{1j}^{(k-1)} ~-~ <1,1>_{j-1}a_{ij}^{(k-1)}}
\end{displaymath}
and let
\begin{displaymath}
d(j) ~\equiv~ <1,1>_{j}^{1/2}.
\end{displaymath}
Finally, let
\begin{displaymath}
\hat{a}_{ij}^{(k)} = 
\left\{
\begin{array}{ll}
0&\mbox{if $i=1$, $1 < j \leq n-k+1$}\\*[.15in]
<1,i>_{n-k+1}/{d(n-k+1)}&
\mbox{if $1 \leq i \leq n-k+1$, $j = 1$}\\*[.15in]
{n(i,j,k)}/{d(j-1) d(j)}&
\mbox{if $1 < i,j \leq n-k+1$}\\*[.15in]
a_{ij}^{(k-1)}&\mbox{otherwise}.
\end{array}
\right.
\end{displaymath}
Then $B = A\Theta = A^{(n)}$. For example, suppose that
\begin{displaymath}
A = A^{(0)} =
\left[
\begin{array}{rrr}
1&1&1\\
0&1&1\\
0&0&1
\end{array}
\right]
\end{displaymath} 
Applying Equation (\ref{eqn:problem}) gives
\begin{displaymath}
A^{(1)} =
\left[
\begin{array}{rrr}
-1/\sqrt{2}&-1/\sqrt{2}\sqrt{3}&2/\sqrt{3}\\
0&-\sqrt{2}/\sqrt{3}&1/\sqrt{3}\\
0&0&\sqrt{3}
\end{array}
\right]
\end{displaymath}
\begin{displaymath}
A^{(2)} =
\left[
\begin{array}{rrr}
1/\sqrt{2}&1/\sqrt{3}&\sqrt{3}/3\sqrt{2}\\
0&\sqrt{3}&0\\
0&2/\sqrt{3}&\sqrt{2}/\sqrt{3}
\end{array}
\right]
\end{displaymath}
\begin{displaymath}
A^{(3)} = 
\left[
\begin{array}{rrr}
\sqrt{3}&0&0\\
2/\sqrt{3}&\sqrt{2}/\sqrt{3}&0\\
1/\sqrt{3}&\sqrt{3}/3\sqrt{2}&1/\sqrt{2}
\end{array}
\right]
\end{displaymath} 
The reader can verify that $B = A^{(3)}$ by using the identity
$BB^T = AA^T$.

\subsubsection{Matrix Transpose}
We are given an $n \times n$ matrix $A = \{ a_{ij} \}$. The problem is
to compute $B = A^{T}$. In Equation (\ref{eqn:problem}), let $m=n$. Let
\begin{displaymath}
\hat{a}_{ij}^{(k)} = 
\left\{
\begin{array}{ll}
a_{ji}^{(k-1)}&\mbox{if $i=1$ and $1 < j \leq n-k+1$}\\*[.1in]
a_{ji}^{(k-1)}&\mbox{if $1 < i \leq n-k+1$ and $j=1$}\\*[.1in]
a_{ij}^{(k-1)}&\mbox{otherwise}.
\end{array}
\right.
\end{displaymath}
Then $B = A^{T} = A^{(n)}$. 

\section{Mapping on Linear Array}
In this section we give a simple approach for solving the recurrence
in Equation (\ref{eqn:problem}) using a unidirectional linear array,
and we illustrate the approach using examples from the previous
section.

\subsection{Machine Model}
The machine model is the linear array depicted in Figure
\ref{fig:array}. The model is an informal idealization of the Warp
computer \cite{annaratone87e}.  Rectangles denote computational cells.
Arrows denote unidirectional communications channels. Cells are
numbered from left to right, starting with unity.

Cells are separately programmable and asynchronous.  Each cell has its
own local code and data storage. Each cell knows its number, which is
contained in local variable {\tt cellid}. Each cell runs at most one
process, which is assigned to the cell at compile time. The
communications channel between two cells is a bounded FIFO buffer.
Data passes through the array from left to right.  Cell $k$ sends a
word of data to cell $k+1$ using a statement of the form
\begin{displaymath}
\mbox{\tt {\bf send}(val)}
\end{displaymath}
which blocks until the buffer between cells $k$ and $k+1$ is nonfull,
and then places the local variable {\tt val} at the end of the buffer.
Similarly, cell $k$ receives a word of data from cell $k-1$ using a
statement of the form
\begin{displaymath}
\mbox{\tt {\bf receive}(val)}
\end{displaymath}
which blocks until the buffer between cells $k$ and $k-1$ is nonempty,
and then retrieves the first item from the buffer and assigns it to
variable {\tt val}.

\subsection{Basic Algorithm}
The mapping of matrices and transformations onto the linear array is
shown in Figure \ref{fig:mapping}. Cell $k$ receives $A^{(k-1)}$ in
Equation (\ref{eqn:problem}) from its left hand neighbor, computes the
transformation
\begin{displaymath}
\hat{A}^{(k)} = T^{(k)}(A^{(k-1)}),
\end{displaymath}
permutes $\hat{A}^{(k)}$, and sends the resulting $A^{(k)}$ in
Equation (\ref{eqn:problem}) to its right hand neighbor.  The scalar
elements of $A^{(k-1)}$ are assumed to arrive at cell $k$ in row-major
order, although in general other orderings (such as column-major or
diagonal-major) can be used.

Figure \ref{fig:basic} describes the behavior of cell $k$.  In this
basic algorithm, the transformation $T^{(k)}$ in Equation
(\ref{eqn:problem}) is simply the null transformation
\begin{displaymath}
\hat{A}^{(k)} = A^{(k-1)}.
\end{displaymath}
By applying less trivial transformations on each cell, or equivalently,
by adding some statements to the basic algorithm in Figure
\ref{fig:basic}, we are able to solve any recurrence in the form of
Equation (\ref{eqn:problem}), including all of the examples from the
previous section.

For example, in Figure \ref{fig:sp}, adding Line (7a) to the basic
algorithm gives a shortest paths algorithm. Notice that we take
advantage of the invariance of the elements in the first row and
column of the input matrix.  In Figure \ref{fig:mi}, adding lines
(2a)--(2c), (5a), (7a), and (12a) to the basic algorithm gives a
matrix inverse algorithm. Finally, in Figure \ref{fig:gr}, adding a
number of lines to the basic algorithm gives an algorithm for matrix
triangularization.

The approach we have described offers a number of practical
qualitative advantages. The first advantage is uniformity of design.
Once we are able to cast a problem in the form of Equation
(\ref{eqn:problem}), the mapping onto the linear array is
straightforward and mechanical. We simply take the basic algorithm in
Figure \ref{fig:basic} and add the code that implements the
transformation $T^{(k)}$. The basic pattern of {\tt \bf send} and {\tt
\bf receive} statements (always a tricky part of a design) does not
change!  The second advantage is simplicity of programming. Each cell
in the linear array executes the same code; we are not forced to
program each cell separately.  The third advantage is that the
approach generalizes easily to large problems on arrays with fewer
than $m$ cells. As a simple example, suppose that we want to compute
the shortest paths for a graph with $n = 350$ vertices using an array
with 10 cells. If we stream the initial adjacency matrix $A$ through
the array, the output of cell 10 is $A^{(10)}$.  If we then stream
$A^{(10)}$ through the array, the output of cell 10 is $A^{(20)}$. And
so on. In fact, this is precisely the technique used by P. S.  Tseng
\cite{pstseng87} at Carnegie Mellon for the implementation of the
shortest paths problem in \cite{annaratone87e}.

Assume that the transformation of each scalar element of $A^{(k-1)}$
can be applied in constant time on cell $k$ (true for all problems
discussed in this paper). Then cell $k$ waits $O(kn)$ time for the
first element of $A^{(k-1)}$ to arrive from the left, and then
transforms and permutes the elements of $A^{(k-1)}$ as they arrive in
$O(n^2)$ time. Thus, using $m$ cells, the result matrix $B = A^{(m)}$
in Equation (\ref{eqn:problem}) is computed in $O(n^2 + mn)$ time.
This points out a major quantitative disadvantage of the approach; it
does not necessarily result in an optimal algorithm. For example, the
parallel algorithm for matrix multiplication in \cite{annaratone87e}
requires only $O(n^{2})$ time on a linear array.  As another
example, the {\em sequential} algorithm for matrix transpose requires
only $O(n^{2})$ time.

\section{Summary and Discussion} 
We have presented two results. First, we demonstrated that a number of
classical problems from linear algebra and graph theory --- in
particular, instances of the algebraic path problem, matrix
multiplication, matrix triangularization, and matrix transpose --- can be
solved using the same recurrence. Second, we gave a simple mapping of
the recurrence onto a linear array.

We have reported on an interesting phenomenon, namely, that a number
of interesting problems can be solved on a linear array in a similar
way.  Unfortunately, we do not yet understand the underlying
mathematical characteristic of these problems that gives rise to this
phenomenon.  We identified these problems more by chance, as we
developed algorithms and programs for a prototype of the Warp
computer, than by any consistent and organized search.  The Warp
prototype supported only unidirectional communication, with the same
code executing on each of its 10 cells. In some cases, the approach
presented in this paper was our only alternative.  Essentially, we
noticed that for many problems, the $kth$ row and column of the matrix
are in some way ``special'' during the $kth$ iteration over the
matrix.  Further, the values of the $kth$ row and column are needed to
compute other elements of the matrix during the $kth$ iteration.  The
intuition for the approach we have presented is that the permutation
of the matrix in Equation (\ref{eqn:problem}) keeps these special rows
and columns in the same position for each cell in the linear array,
and the restriction on transformation $T^{(k)}$ in Equation
(\ref{eqn:problem}) guarantees that all of the needed matrix elements
have been received by a cell when it needs to pass some result on to
the next cell.

The fact that a single column and row of the matrix are used to
compute the other elements of the matrix is reminiscent of the rank
one modification from linear algebra. Perhaps the underlying
mathematical characteristic of the problems discussed in this paper
can ultimately be explained in these terms. 

\section*{Acknowledgements}
Thanks to H. T. Kung for providing access to the Warp machine at
Carnegie Mellon. Thanks also to the anonymous referees, and to Chuck
Fiduccia, Bob Mattheyses, and P. S. Tseng for helpful discussions.

\bibliography{/usr/droh/bib/defs,/usr/droh/bib/refs,/usr/droh/bib/me}

\begin{figure}[p]
\setlength{\unitlength}{0.2in}
\begin{center}
\begin{picture}(15,30)(1,4)
%vertices
\multiput(4,13)(2,0){4}{\makebox(0,0){$\circ$}}
\multiput(4,11)(2,0){4}{\makebox(0,0){$\circ$}}
\multiput(4,9)(2,0){4}{\makebox(0,0){$\circ$}}
\multiput(4,7)(2,0){4}{\makebox(0,0){$\circ$}}
\multiput(7,10)(2,0){4}{\makebox(0,0){$\bullet$}}
\multiput(7,8)(2,0){4}{\makebox(0,0){$\bullet$}}
\multiput(7,6)(2,0){4}{\makebox(0,0){$\bullet$}}
\multiput(7,4)(2,0){4}{\makebox(0,0){$\bullet$}}
%circle edges
\multiput(4.1,13)(2,0){3}{\vector(1,0){1.8}}
\multiput(4.1,11)(2,0){3}{\vector(1,0){1.8}}
\multiput(4.1,9)(2,0){3}{\vector(1,0){1.8}}
\multiput(4.1,7)(2,0){3}{\vector(1,0){1.8}}
\multiput(4,12.9)(2,0){4}{\vector(0,-1){1.8}}
\multiput(4,10.9)(2,0){4}{\vector(0,-1){1.8}}
\multiput(4,8.9)(2,0){4}{\vector(0,-1){1.8}}
%bullet edges
\multiput(7.1,10)(2,0){3}{\vector(1,0){1.8}}
\multiput(7.1,8)(2,0){3}{\vector(1,0){1.8}}
\multiput(7.1,6)(2,0){3}{\vector(1,0){1.8}}
\multiput(7.1,4)(2,0){3}{\vector(1,0){1.8}}
\multiput(7,9.9)(2,0){4}{\vector(0,-1){1.8}}
\multiput(7,7.9)(2,0){4}{\vector(0,-1){1.8}}
\multiput(7,5.9)(2,0){4}{\vector(0,-1){1.8}}
%diagonal edges
\multiput(6.1,10.9)(2,0){3}{\vector(1,-1){0.8}}
\multiput(6.1,8.9)(2,0){3}{\vector(1,-1){0.8}}
\multiput(6.1,6.9)(2,0){3}{\vector(1,-1){0.8}}
%labels
\put(2.5,10){\makebox(0,0){$A^{(k-1)}$}}
\put(14,7){\makebox(0,0){$A^{(k)}$}}
\end{picture}
\end{center}
\caption{Visual interpretation of set $S_{ij}$ for $n=4$}
\label{fig:pgraph}
\end{figure}

\begin{figure}[p]
\setlength{\unitlength}{0.2in}
\begin{center}
\begin{picture}(6,30)(0,1)
%vertices
\put(2,5){\makebox(0,0){$\bullet$}}
\put(2,2){\makebox(0,0){$\bullet$}}
\put(5,5){\makebox(0,0){$\bullet$}}
\put(5,2){\makebox(0,0){$\bullet$}}
\put(2,5.5){\makebox(0,0){$v_1$}}
\put(2,1.5){\makebox(0,0){$v_2$}}
\put(5,5.5){\makebox(0,0){$v_3$}}
\put(5,1.5){\makebox(0,0){$v_4$}}
%edges
\put(2,5){\line(1,0){3}}
\put(2,2){\line(1,0){3}}
\put(2,5){\line(0,-1){3}}
\put(5,5){\line(0,-1){3}}
\put(1.3,3.5){$2$}
\put(3.5,5.3){$1$}
\put(5.3,3.5){$7$}
\put(3.5,1.3){$3$}
\end{picture}
\end{center}
\caption{Graph for algebraic path problem examples, $n=4$}
\label{fig:exgraph}
\end{figure}

\begin{figure}[p]
\setlength{\unitlength}{0.2in}
\begin{center}
\small
\begin{picture}(20,30)(0,0)
%cells
\put(3,1){\framebox(2,2){1}}
\put(7,1){\framebox(2,2){2}}
\put(12,2){\makebox(0,0){$\cdots$}}
\put(15,1){\framebox(2,2){m}}
%channels
\multiput(1,2)(4,0){5}{\vector(1,0){2}}
\end{picture}
\end{center}
\caption{Linear array machine model}
\label{fig:array}
\end{figure}

\begin{figure}[p]
\setlength{\unitlength}{0.2in}
\begin{center}
\begin{picture}(20,30)(0,0)
%cells
\put(3,1){\framebox(2,2){$T^{(1)}$}}
\put(7,1){\framebox(2,2){$T^{(2)}$}}
\put(12,2){\makebox(0,0){$\cdots$}}
\put(15,1){\framebox(2,2){$T^{(m)}$}}
%channels
\multiput(1,2)(4,0){5}{\vector(1,0){2}}
%data
\put(2,2.2){\makebox(0,0)[b]{$A^{(0)}$}}
\put(6,2.2){\makebox(0,0)[b]{$A^{(1)}$}}
\put(10,2.2){\makebox(0,0)[b]{$A^{(2)}$}}
\put(13.8,2.2){\makebox(0,0)[b]{$A^{(m-1)}$}}
\put(18,2.2){\makebox(0,0)[b]{$A^{(m)}$}}
\end{picture}
\end{center}
\caption{Mapping of the recurrence onto a linear array}
\label{fig:mapping}
\end{figure}

\begin{figure}[p]
\baselineskip=5pt
\vspace{4in}
\begin{center}
\begin{tabbing}
\tt
12345123451234512345\=12345\=12\=12\=12\=12\=12\=12\=12\=12\=12\=12\= \kill
\+ \kill
1.\>{\bf for} j := 1 {\bf to} n {\bf do}\\
2.\>\>{\bf receive}(row[j]);\\
3.\>{\bf endfor}\\
4.\>{\bf for} i := 2 {\bf to} n {\bf do}\\
5.\>\>{\bf receive}(col);\\
6.\>\>{\bf for} j := 2 {\bf to} n {\bf do}\\
7.\>\>\>{\bf receive}(val);\\
8.\>\>\>{\bf send}(val);\\
9.\>\>{\bf endfor}\\
10.\>\>{\bf send}(col);\\
11.\>{\bf endfor}\\
12.\>{\bf for} j := 2 {\bf to} n {\bf do}\\
13.\>\>{\bf send}(row[j]);\\
14.\>{\bf endfor}\\
15.\>{\bf send}(row[1]);\\
\end{tabbing}					
\caption{Algorithm for cell $k$ (basic algorithm)}
\label{fig:basic}
\end{center}
\end{figure}					

\begin{figure}[p]
\vspace{4in}
\begin{center}
\begin{tabbing}
\tt
12345123451234512345\=12345\=12\=12\=12\=12\=12\=12\=12\=12\=12\=12\= \kill
\+\kill
1.\>\>{\bf for} j := 1 {\bf to} n {\bf do}\\
2.\>\>\>{\bf receive}(row[j]);\\
3.\>\>{\bf endfor}\\
4.\>\>{\bf for} i := 2 {\bf to} n {\bf do}\\
5.\>\>\>{\bf receive}(col);\\
6.\>\>\>{\bf for} j := 2 {\bf to} n {\bf do}\\
7.\>\>\>\>{\bf receive}(val);\\
7a.\>\>\>\>val := {\bf min}(val, col + row[j]);\\
8.\>\>\>\>{\bf send}(val);\\
9.\>\>\>{\bf endfor}\\
10.\>\>\>{\bf send}(col);\\
11.\>\>{\bf endfor}\\
12.\>\>{\bf for} j := 2 {\bf to} n {\bf do}\\
13.\>\>\>{\bf send}(row[j]);\\
14.\>\>{\bf endfor}\\
15.\>\>{\bf send}(row[1]);
\end{tabbing}					
\caption{Algorithm for cell $k$ (shortest paths)}
\label{fig:sp}
\end{center}
\end{figure}

\begin{figure}[p]
\vspace{4in}
\begin{center}
\begin{tabbing}
\tt
12345123451234512345\=12345\=12\=12\=12\=12\=12\=12\=12\=12\=12\=12\= \kill
\+\kill
1.\>\>{\bf for} j := 1 {\bf to} n {\bf do}\\
2.\>\>\>{\bf receive}(row[j]);\\
2a.\>\>\>{\bf if} i = 1 {\bf then}\\ 
2b.\>\>\>\>row[1] := 1/(1 - row[1]);\\
2c.\>\>\>{\bf endif}\\
3.\>\>{\bf endfor}\\
4.\>\>{\bf for} i := 2 {\bf to} n {\bf do}\\
5.\>\>\>{\bf receive}(col);\\
5a.\>\>\>col := col * row[1];\\
6.\>\>\>{\bf for} j := 2 {\bf to} n {\bf do}\\
7.\>\>\>\>{\bf receive}(val);\\
7a.\>\>\>\>val := val + row[j]*col;\\
8.\>\>\>\>{\bf send}(val);\\
9.\>\>\>{\bf endfor}\\
10.\>\>\>{\bf send}(col);\\
11.\>\>{\bf endfor}\\
12.\>\>{\bf for} j := 2 {\bf to} n {\bf do}\\
12a.\>\>\>row[j] := row[j]*row[1];\\
13.\>\>\>{\bf send}(row[j]);\\
14.\>\>{\bf endfor}\\
15.\>\>{\bf send}(row[1]);
\end{tabbing}					
\caption{Algorithm for cell $k$ (matrix inverse)}
\label{fig:mi}
\end{center}
\end{figure}					

\begin{figure}[p]
\begin{center}
\begin{tabbing}
\tt
12345678901234\=12345\=12\=12\=12\=12\=12\=12\=12\=12\=12\=12\= \kill
\+\kill
\>\>{\bf initially:}\\
\>\>ip11[0] := 0;\\ 
\>\>last := n-cellid+1;\\*[\baselineskip]
1.\>\>{\bf for} j := 1 {\bf to} n {\bf do}\\
2.\>\>\>{\bf receive}(row[j]);\\
2a.\>\>\>{\bf if} j $\leq$ last {\bf then}\\
2b.\>\>\>\>ip11[j] := ip11[j-1] + row[j]*row[j];\\
2c.\>\>\>{\bf endif}\\
3.\>\>{\bf endfor}\\
4.\>\>{\bf for} i := 2 {\bf to} n {\bf do}\\
5.\>\>\>{\bf receive}(col);\\
5a.\>\>\>ip1i[1] := col*row[1];\\
6.\>\>\>{\bf for} j := 2 {\bf to} n {\bf do}\\
7.\>\>\>\>{\bf receive}(val);\\
7a.\>\>\>\>{\bf if} j $\leq$ last {\bf then}\\
7b.\>\>\>\>\>ip1i[j] := ip1i[j-1] + val*row[j];\\
7c.\>\>\>\>\>numer := ip1i[j-1]*col - ip11[j-1]*val;\\
7d.\>\>\>\>\>denom := {\bf sqrt}(ip11[j-1])*{\bf sqrt}(ip11[j]);\\
7e.\>\>\>\>\>val := numer/denom;\\
7f.\>\>\>\>{\bf endif}\\
8.\>\>\>\>{\bf send}(val);\\
9.\>\>\>{\bf endfor}\\
9a.\>\>\>{\bf if} i $\leq$ last {\bf then}\\
9b.\>\>\>\>col := ip1i[last]/{\bf sqrt}(ip11[last]);\\
9c.\>\>\>{\bf endif}\\
10.\>\>\>{\bf send}(col);\\
11.\>\>{\bf endfor}\\
12.\>\>{\bf for} j := 2 {\bf to} n {\bf do}\\
12a.\>\>{\bf if} j $\leq$ last {\bf then}\\
12b.\>\>\>row[j] := 0;\\
12c.\>\>{\bf endif}\\
13.\>\>\>{\bf send}(row[j]);\\
14.\>\>{\bf endfor}\\
14a.\>\>row[1] := {\bf sqrt}(ip11[last]);\\
15.\>\>{\bf send}(row[1]);
\end{tabbing}					
\caption{Algorithm for cell $k$ (matrix triangularization)}
\label{fig:gr}
\end{center}
\end{figure}					

\newpage
\listoffigures
\newpage
\section*{Mailing Address and Phone Number}
Please return galley proofs to\\
\\
David R. O'Hallaron\\
School of Computer Science\\
Carnegie Mellon University\\
Pittsburgh, PA, 15213\\
\\
phone:(412)268-8199

\newpage

\section*{Biography}

David R. O'Hallaron received the Ph.D. degree in Computer Science from
the University of Virginia, Charlottesville, in 1986.

From 1986 to 1989 he was with the GE Research and Development Center,
Schenectady, NY, as a Computer Scientist. In 1989, he joined the
faculty in the School of Computer Science, Carnegie Mellon University,
Pittsburgh, PA, as a Systems Scientist. His main research interests
are in the areas of programming models, algorithms, and applications
for parallel machines.

\end{document}
