\makeatletter
\long\def\unmarkedfootnote#1{{\long\def\@makefntext##1{##1}\footnotetext{#1}}}
\makeatother
\newcommand{\comment}[1]{}

\documentstyle[epsf]{cmu-art}
\begin{document}
\bibliographystyle{acm}

\section{Example applications}
In this section, we apply the framework developed in the previous
section to three realistic task and data parallel applications: fast
Fourier transform~\cite{subhlok93}, narrowband tracking
radar~\cite{shaw92}, and multibaseline stereo imaging~\cite{webb93}.
We briefly describe each application, show how it is implemented as an
Fx program, and use the framework to reason about mapping task
parallel programs on a multicomputer. Each application has unique
properties that excercises different parts of the mapping framework.
We have found these applications to be invaluable; interested readers
can obtain complete Fortran 77 implementations from the authors.

\subsection{Fast Fourier transform}
The fast Fourier transform (FFT) is of fundamental importance in
signal and image processing.  The FFT is interesting to study because
its task parallel structure is similar to many larger signal
processing applications, and because it is a satisfying combination of
the easy and the challenging: easy because there is perfect data
parallelism within each task; challenging because it must scan both
the rows {\em and} the columns of the input array, and therefore
requires an efficient redistribution of data to achieve good
performance.

The FFT is an efficient procedure for computing the discrete Fourier
transform (DFT). A 1D~DFT is a complex matrix-vector product
\begin{displaymath}
y = F_n x
\end{displaymath}
where $x = [x_0 , \ldots , x_{n-1}]$ and $y = [y_0 ,\ldots, y_{n-1}]$
are complex vectors, and $F_n = (f_{pq})$ is an $n \times n$ matrix
such that
\begin{displaymath}
f_{pq} = \omega_n^{pq}, ~~0 \leq p,q \leq n-1, 
\end{displaymath}
where
\begin{displaymath}
\omega_n = \cos(2\pi/n) - i \sin(2\pi/n) = e ^ {-2 \pi i /n }
\end{displaymath}
and $i = \sqrt{-1}$.  An FFT algorithm exploits structure in $F_n$ to
compute the DFT in $O(n \log n)$ floating point operations.  See
\cite{vanloan92} for an excellent description of the numerous FFT
algorithms that have been developed over the past 40 years.

If $n=n_1 n_2$, then a 1D~FFT can be computed using a collection of
smaller independent 1D~FFTs and scalings
\cite{bailey90,gentleman66,vanloan92}.
Starting with a view of the input vector $x$ as an $n_1 \times n_2$
matrix $A$, stored in column major order, perform $n_1$ $n_2$-point
FFTs on the rows of $A$, scale each element $a_{pq}$ of the resulting
matrix by a factor of $\omega_n^pq$, and then perform $n_2$
$n_1$-point FFTs on the columns of $A$.  Obtain the result vector $y$
by concatenating the rows of the resulting matrix. 
Higher dimensional FFTs are computed in a similar fashion.  Given an
$n_1 \times n_2$ input matrix $A$, perform $n_1$ $n_2$-point 1D~FFTs
on the rows, followed by $n_2$ $n_1$-point FFTs on the columns of the
resulting matrix.  
Figure \ref{fig:fft}(b) summarizes the procedure.
\begin{figure}[htb]
%\epsfxsize=4in
\centerline{\epsfbox{fft.eps}}
\caption{2D FFT benchmark, $n_1 = n_2 = 128$}
\label{fig:fft}
\end{figure}

The issues involved in mapping 1D and 2D FFTs are similar. Here we will
focus on  the issues involved in mapping a $128 \times 128$ 2D~FFT on
64 processors. In Fx, an FFT program that operates on a stream of $m$
arrays has the following form:
\tt
\begin{tabbing}
tttttttt\=ttt\=ttt\= ttt\=ttt\= \kill
C\$ \>{\bf begin parallel} \\
\>     do i = 1,m \\
\>\>       call rowffts(A)\\
C\$\>\>     {\bf output:} A\\
\>\>       call colffts(A)\\
C\$\>\>     {\bf input:} A\\
\>     enddo\\
C\$\>   {\bf end parallel}
\end{tabbing}
\rm
\normalsize
The program consists of a parallel section, with calls to two
task-subroutines inside a loop that iterates $m$ times.  Subroutine
{\tt rowffts} is a data parallel procedure that produces a complex
$128 \times 128$ array $A$, distributed by rows, and applies a 1D~FFT
to each row.  Subroutine {\tt rowffts} is a data parallel procedure
that redistributes $A$ by columns, and applies a 1D~FFT to each
column.  The source of data parallelism within each task-procedure is
a parallel loop statement. 

\subsection*{Mapping the FFT program}
Following the framework described in Section \ref{sec:framework}, we
first decide how to partition the two task-subroutines, {\tt rowffts}
and {\tt colffts}. There are two possibilities: assign the
task-subroutines to the same module or to two different modules.
According to the framework, we evaluate the effectiveness of the two
stategies by applying the following tests:

First we consider the scalability of the individual tasks.  The
data-parallelism inside each task-subroutine is in the form of a
parallel loop with no dependences between loop iterations.  Thus each
task-subroutine scales linearly with the number of processors.  This
program property is captured by the curve labeled ``Perfectly scalable
computation'' in Figure~\ref{fig:xxx}.  Because each task scales well,
there is no significant benefit to either partitioning strategy,
according to this test. Notice also that this implies no significant
benefit to replication either.

Next we consider the memory requirement and decoupling tests.  But
since replication is not beneficial, these tests do not help us reach
a conclusion either.

Finally, we consider the costs of the inter-task communication
associated with the two partitioning strategies. When the 2D~FFT
program transfers an $n \times n$ array from {\tt rowffts} to {\tt
colffts}, it must redistribute the array from a row-wise distribution
to a column-wise distribution.  For the iWarp version of the Fx
system, and for roughly $n \leq 256$, this redistribution is cheaper
(for known and predictable technical reasons) if it occurs between two
modules, each running on 32 nodes, than if it occurs on one module
running on 64 nodes. So this test suggests partitioning into two
modules.

Next, according to the framework, we consider the replication of
modules.  Since both task-procedures in the 2D~FFT program are
perfectly scalable, there is little or no benefit to replication.  

Based on this analysis, we decide to partition the task-subroutines
into two unreplicated modules, as shown in Figure~\ref{fig:ffttask}.  The
final mapping, after processor allocation, consists of two modules,
each running on 32 nodes.  Empirical results confirm that this mapping
is indeed more efficient for smaller data sets~\cite{subhlok93}.

\begin{figure}[htb]
%\epsfxsize=4in
\centerline{\epsfbox{ffttask.eps}}
\caption{Mapping the FFT}
\label{fig:ffttask}
\end{figure}

\subsection{Narrowband tracking radar}
The narrowband tracking radar benchmark was developed by researchers
at MIT Lincoln Labs to measure the effectiveness of various
multicomputers for their radar applications~\cite{shaw92}. It is a
particularly interesting benchmark for studying task parallelism
because of hard real-time requirements, and because the size of the
input data set is limited by physical properties of the radar unit.
The amount of available low-level data parallelism is limited, so
additional parallelism must come from higher-level task parallelism.

The radar algorithm is shown in Figure~\ref{fig:radar}.  The program
inputs data from a single sensor along $c=4$ independent {\em
channels}.  Every 5 milliseconds, for each channel, the program
receives $d=512$ complex vectors of length $r=10$, one after the other
in the form of an $r \times d$ complex matrix $A$ (assuming the column
major ordering of Fortran).  At a high-level, each input matrix $A$ is
processed in the following way: (1) {\em Corner turn} the $r \times d$
input matrix to form a $d \times r$ matrix.  (2) Perform $r$
independent $d$-point FFTs.  (3) Convert the resulting complex $d
\times r$ matrix to a real $w \times r$ submatrix, $w=40$, by
replacing each element $a + ib$ in the $w \times r$ submatrix with its
scaled magnitude $\sqrt{a^2 + b^2}/d$.  (4) Threshold each element
$a_{jk}$ of the submatrix using a cutoff that is a function of
$a_{jk}$ and the sum of the submatrix elements.

\begin{figure}[htb]
%\epsfxsize=4in
\centerline{\epsfbox{radar.eps}}
\caption{Radar benchmark, $c=4$, $r=10$, $d=512$}
\label{fig:radar}
\end{figure}

The Fx version of the radar program operating on a stream of $m$ input
data sets has the following form:
\tt
\begin{tabbing}
tttttttt\=ttt\=ttt\= ttt\=ttt\= \kill
C\$ \>{\bf begin parallel} \\
\>     do i = 1,m \\
\>\>       call datagen(A)\\
C\$\>\>     {\bf output:} A\\
\>\>       call compute(A,B)\\
C\$\>\>     {\bf input:} A\\
C\$\>\>     {\bf output:} B\\
\>     enddo\\
C\$\>   {\bf end parallel}
\end{tabbing}
\rm
\normalsize
Task-subroutine {\tt datagen} runs on one processor and acquires the
data from all 4 channels. Task-subroutine {\tt compute} is a data
parallel routine that corner turns the input array, performs
independent sequential FFTs on the columns, scales the result,
performs thresholding, and returns the result.  The data-parallelism
inside the task-subroutine is in the form of a parallel loop where
each loop iteration operates on a single column.  

\subsection*{Mapping the radar program}
We must first decide how to partition the two task-subroutines, {\tt
datagen} and {\tt compute}. As with the FFT example, there are two
possibilities: assign the task-subroutines to the same module or to
two different modules. Again, we evaluate the effectiveness of the two
strategies by applying four tests:

First, we consider the scalability of the individual task-subroutines.
{\tt Datagen} acquires input from the sensor at a constant rate, so it
does not scale well.  Because there are only 10 columns in its input
array, {\tt compute} scales perfectly up to 10 processors, then drops
off rapidly.  This type of behavior is captured by the curve labeled
``Scaleable up to p processors'' in Figure~\ref{fig:xxx}.  Because
{\tt datagen} scales so poorly, there is little benefit to
partitioning the tasks onto two modules, according to this test.

The memory requirement test is also inconclusive because the data sets
are small.  However the decoupling test provides some insight.  Since
each instance of {\tt compute} can use at most 10 processors
effectively, replication might provide a way to utilize more
processors. Since {\tt datagen} runs efficiently on only one
processor, this test suggests partitioning the task-subroutines into
two modules so that {\tt compute} can be replicated in a later step.

Finally we consider the relative costs of inter-task communication for
the two partitioning strategies.  As we noted for the FFT program,
inter-task communication of small data sets is more efficient between
separate modules than within the same module on our system.  So this
test also suggests partitioning the task-subroutines into two modules.

Next we consider the replication of modules.  We have already
determined that replication of the {\tt compute} module could be
beneficial. So let us consider the possible constraints on replication
identified in the framework.  There are no true dependences between
different {\tt compute} tasks running on separate modules, so
replication is allowed. Next we consider constraints on the amount of
replication. The data sets are small and could fit on one processor,
so memory requirements are not an issue.  However, on our system,
limited communication resources restrict the replication to four
instances of the {\tt compute} module. Another possible constraint is
the tolerance of the application for higher latencies.  For the radar
example, latencies introduced by 4-way replication are probably
tolerable.

The decision then is to partition the task-subroutines into two
modules, and to replicate the {\tt compute} module four times. After
processor allocation the final mapping consists of a {\tt datagen}
module running on one processor, and 4 {\tt compute} modules, each
running on ten processors.  This is shown in Figure~\ref{fig:radartask}.

\begin{figure}[tbh]
%\epsfxsize=4in
\centerline{\epsfbox{radartask.eps}}
\caption{Radar task graph and mapping}
\label{fig:radartask}
\end{figure}

<paragraph giving performance of data-parallel, replicated-2,
and replicated-4 cases>

\subsection{Multibaseline stereo}
This program, which was developed at Carnegie Mellon, uses the
information from multiple video cameras to display depth information
in real time~\cite{webb93}. It is an interesting program because,
unlike the FFT and radar examples, there are data-parallel tasks that
contain a significant amount of intra-task communication.

The basic stereo algorithm is shown in Figure~\ref{fig:stereo}.  
\begin{figure}[htb]
%\epsfxsize=4in
\centerline{\epsfbox{stereo.eps}}
\caption{Multibaseline stereo benchmark, $m=480$, $n=512$}
\label{fig:stereo}
\end{figure}
Input consists of three $m \times n$ images acquired from three
horizontally aligned, equally spaced cameras. One image is the {\em
reference image}, the other two are {\em match images}. For each of 16
disparities, $d = 1,\ldots,16$, the first match image is shifted by
$d$ pixels, the second image is shifted by $2d$ pixels. A {\em
difference image} is formed by computing the sum of squared
differences between the pixels of the reference image and the shifted
match images.  Next, an {\em error image} is formed by replacing each
pixel in the difference image with the sum of the pixels in a
surrounding $13 \times 13$ window. A {\em disparity image} is then
formed by finding, for each pixel, the disparity that minimizes error.
The depth of each pixel is then be displayed as a simple function of
its disparity.

The Fx version of the stereo  program operating on a stream of $s$ input
data sets has the following form:
\tt
\begin{tabbing}
tttttttt\=ttt\=ttt\= ttt\=ttt\= \kill
C\$ \>{\bf begin parallel} \\
\>     do i = 1,s \\
\>\>       call datagen(R,M1,M2)\\
C\$\>\>    {\bf output:} R,M1,M2\\
\>\>       do d = 1,16\\
\>\>\>		call diff(R,M1,M2,DIFF)\\
C\$\>\>\>       {\bf input:} R,M1,M2\\
C\$\>\>\>       {\bf output:} DIFF\\
\>\>\>		call error(DIFF,ERR(1,1,d)\\
C\$\>\>\>       {\bf input:} DIFF\\
C\$\>\>\>       {\bf output:} ERR(1,1,d)\\
\>\>       enddo\\
\>\>       call reduce(ERR,DISP)\\
C\$\>\>    {\bf input:} ERR\\
C\$\>\>    {\bf output:} DISP\\
\>     enddo\\
C\$\>   {\bf end parallel}
\end{tabbing}
\rm
\normalsize
Task-subroutine {\tt datagen} runs on one processor and acquires the
images from the cameras. Each of the 16 instances of the {\tt diff}
task-subroutine is a perfectly data-parallel routine that converts the
input images to a difference image. Each instance of the {\tt error}
task-subroutine is a data-parallel routine that sums over a window of
pixels in the difference image to produce an error measure for each
pixel. Each image is distributed by rows within each task, so a
processor needs to exchange rows with it neighbors before the error
image can be produced. The output of the {\tt error} tasks is passed
to {\tt reduce}, which applies a min-reduction to produce the
disparity image, and then displays the corresponding depth
information.  

\subsection*{Mapping the stereo program}
As before, we evaluate different partitioning strategies using the
four tests in the framework. Although the stereo program is much
richer than the FFT and radar programs, the line of reasoning in
determining a mapping is similar.

First, we consider the scalability of the individual task-subroutines.
{\tt Datagen} acquires input from the cameras at a constant rate and
thus does not scale well. {\tt Diff} scales perfectly because there is
no intra-task communication. Because of significant intra-task
communication, {\tt error} runs less efficiently as more processors
are added, as in the curve labeled ``Computation with communication''
in Figure~\ref{fig:xxx}. {\tt Reduce} behaves similarly. The
conclusion is that we will try to assign each {\tt error} and {\tt
reduce} task to a separate module.

Since more modules implies fewer processors, and thus less memory per
module, and since the input images are relatively large, the memory
requirement test places a constraint on the number of {\tt error}
modules that can be created. The bound on the number of modules, four
for our system, is determined at processor allocation time.

The decoupling test is not relevant for the stereo program.
Replication is not needed because of the large number of tasks
executed by each iteration of the parallel section.

Finally we consider the relative costs of inter-task communication for
different partitioning strategies. The interactions here
are quite subtle and not fully understood. However, one case is
quite clear: assigning the {\tt diff} and {\tt error} task-subroutines
for the same disparity to the same module results in an inter-task
communication cost of zero. 

Replication is not needed for the stereo program, so the final mapping
partition is to assign {\tt datagen} to a module, to assign {\tt
diff}/{\tt error} pairs to some number of modules (4 in this particular
instance), and to assign {\tt
reduce} to a module.  This is shown in
Figure~\ref{fig:stereotask}.

\begin{figure}[htb]
%\epsfxsize=4in
\centerline{\epsfbox{stereotask.eps}}
\caption{Stereo task graph and mapping}
\label{fig:stereotask}
\end{figure}


<paragraph summarizing performance of data parallel versus task parallel>

\bibliography{/afs/cs/project/iwarp/member/droh/bib/compiler,/afs/cs/project/iwarp/member/droh/bib/defs,/afs/cs/project/iwarp/member/droh/bib/me,/afs/cs/project/iwarp/member/droh/bib/refs}
\end{document}


