\newcommand{\comment}[1]{}
\input epsf
\input{psadobe}
\makeatletter
\long\def\unmarkedfootnote#1{{\long\def\@makefntext##1{##1}\footnotetext{#1}}}
\makeatother
%\documentstyle[fleqn,prog]{art-nocolon}
\documentstyle[prog,times]{art-nocolon}
\include{psfig}
%\include{caption}


\textheight = 9in
\textwidth = 7in              % assume output has this much centered
%\voffset = -.5 in
%\hoffset = -1.0in
%
%
\oddsidemargin =  -.25in
\evensidemargin = -.25in
\topmargin = -.5 in

\renewcommand{\floatpagefraction}{1}
\renewcommand{\topfraction}{1}
\renewcommand{\bottomfraction}{1}
\renewcommand{\textfraction}{0}
%\renewcommand{\baselinestretch}{1}

\newtheorem{definition}{Definition}
\newtheorem{theorem}{Theorem}
\newtheorem{lemma}{Lemma}
\newcommand{\ignore}[1]{}

%\textfloatsep = 0.2in
%\intextsep = 0.2in
%\newcommand{\Tau}{{\large\tau\normalsize\mit}}
\newcommand{\Tau}{\tau}
\newcommand{\Mu}{\mu}
\newcommand{\Chi}{{\cal X}}

\begin{document}
\bibliographystyle{acm}
\title{Communication and Memory Requirements as the Basis for \\
       Mapping Task and Data Parallel Programs }
\author{Jaspal Subhlok, David R. O'Hallaron,  Thomas Gross, Peter 
        A. Dinda and Jon Webb \\ 
{\tt \{jass,droh,trg,pdinda,webb+\}@cs.cmu.edu}\\
School of Computer Science\\
Carnegie Mellon University\\
Pittsburgh, PA 15213 
\vspace*{0.2 in} \\
{\em Corresponding Author:} Jaspal Subhlok\\
               Phone: (412)-268-7893\\
               FAX: (412)-681-5739 \\
               Email: jass@cs.cmu.edu\\
}

\thispagestyle{empty}
\maketitle
\abstract{
\normalsize
For a wide variety of applications, both task and data
parallelism must be exploited to achieve the best possible
performance on a multicomputer. Recent research has
underlined the importance of exploiting task and data 
parallelism in a single compiler framework, and such a
compiler can map a single source program in many different
ways onto a parallel machine. There are 
several complex tradeoffs between task and data parallelism,
depending on the characteristics of the program
to be executed, most significantly the memory and communication
requirements, 
and the performance parameters of the target parallel machine. 
In this paper, we isolate and
examine the specific characteristics of  programs
that determine the  performance for different mappings 
on a parallel machine, and present a framework for
obtaining a good mapping. The framework is applicable to applications
that process a stream of input, and whose computation
structure is fairly static and predictable.
We describe three applications 
that were developed with our compiler: fast Fourier transforms,
narrowband tracking radar and multibaseline stereo, examine
the tradeoffs between various mappings for them, and show
how the framework was used to obtain efficient mappings.
We believe that application and tool developers interested
in combining task and data parallelism will find the ideas and
the results in this paper useful.
}
\normalsize
\newpage
\section{Introduction}\label{intro}
Compilation and optimization of programs for private memory parallel
computers has been a very active area of research for several years.
Several parallelizing compilers have been developed for data parallel
programs, including Fortran~D~\cite{TsHK93} and Vienna
Fortran~\cite{ChMZ92}, and for task parallel
programs~\cite{LaRi91,Prin91}.  Recent research shows that a large
class of applications contain task and data parallelism~\cite{Fox91,agrawal93}
and it is important to exploit them in a single compiler
framework~\cite{ChRe93,FoCh92,SSOG93}.
In this paper, we address how program characteristics influence
the tradeoffs in
task and data parallel programs, and present a framework
for finding good parallel mappings.

There are several reasons to support task parallelism, in addition
to data parallelism, in a parallelizing compiler.
For many applications, particularly in the areas of
digital signal processing and image processing,
the problem sizes are relatively small and fixed,
and not enough data parallelism is available
to effectively use the large number of processors
of a massively parallel machine. 
Even when adequate data parallelism is available,
a data parallel mapping may not be the most efficient
way to execute a program, if the program does not scale well.
A complete application often consists of a set of procedures with
different computation and communication requirements, and 
different sets of processors may have to be assigned to different
procedures for the best possible performance.

Using task and data parallelism together,  it is possible to 
map an application in a variety of ways onto a parallel machine.
Consider an application with three coarse grain, pipelined,
data parallel computation stages, with each execution of
a computation stage corresponding to a task.
 A set of mappings for such an application is shown 
in Figure~\ref{fig:forms}. Figure~\ref{fig:forms}(a) shows a
pure data parallel mapping, where 
all processors participate in all computation stages.
Figure~\ref{fig:forms}(b) shows a pure task parallel mapping, where
a subset of processors is dedicated to each computation stage.
It may be possible to have
multiple copies of the data parallel mapping executing on different 
sets of processors, as shown in  Figure~\ref{fig:forms}(c). Finally,
a mix  of task and data parallelism, with replication
is shown in Figure~\ref{fig:forms}(d).

\begin{figure}[p]
\centerline{\psadobenomag{forms.eps}}
%\epsfbox{forms.eps}
\caption{Combinations of Task and Data parallel mappings}
\label{fig:forms}
\end{figure}

The fundamental question addressed in this paper is: "Which of the
many mappings that are feasible when task and data parallelism are
exploited together is optimal, and how can it be determined?"

We examine  three critical measurable characteristics
of task and data parallel programs; scalability, memory requirement
and inter-task communication cost. We build a framework to
methodically choose a task and data parallel mapping, based
on the program characteristics. The central idea is to find
the best use of available memory and communication resources
to minimize global inter-task and intra-task communication overheads.
We show how the framework was
used to map three applications in signal and image processing;
fast Fourier transforms,
narrowband tracking radar and multibaseline stereo imaging;
and present performance results for them.

\section{Programming and compiling task parallelism}\label{tasks}
Fx is a compiler system that supports task and data
parallelism\cite{SSOG93}.  The base language is Fortran~77 augmented
with Fortran~90 array syntax, and data layout statements based on
Fortran~D\cite{TsHK93} and High Performance Fortran~\cite{HPFF93}.
Data parallelism
is expressed using array syntax and parallel loops.  The main concepts
in the Fx compiler are language independent, and Fortran was chosen
for convenience and user acceptance.

Task parallelism is expressed in special code regions called 
{\em parallel sections}. The body of a parallel section is
limited to calls to subroutines called {\em task-subroutines},
with each execution instance representing a parallel task,
and loops to generate multiple instances.
Each task-subroutine
call is followed by {\tt input} and {\tt output} directives, which
define the interface of the task-subroutine to the calling routine,
that is, they list the  variables in the calling routine
that are accessed and modified  in the called task-subroutine. The
entries in the input and output lists can be scalars, array slices, 
or whole arrays. The task-subroutines may have data parallelism inside
them. A parallel section of an example program is shown in
Figure~\ref{fig:tasks}.

\begin{figure}[p]
\hspace{.1in}
\begin{minipage}[t]{3.0in}
\tt
\begin{tabbing}
ttttttt\=ttt\=ttt\= ttt\=ttt\= \kill
C\$ \>{\bf begin parallel} \\
\>     do i = 1,10 \\
\>\>     call src(A,B)\\
C\$\>\>  {\bf output:} A,B\\
\>\>     call p1(A)\\
C\$\>\>  {\bf input:} A\\
C\$\>\>  {\bf output:} A\\
\>\>     call p2(B)\\
C\$\>\>  {\bf input:} B\\
C\$\>\>  {\bf output:} B\\
\>\>     call sink(A,B)\\
C\$\>\>  {\bf input:} A,B\\
\>     enddo\\
C\$\>   {\bf end parallel}
\end{tabbing}
\rm
\normalsize
\end{minipage} \
\begin{minipage}[t]{4.2in}
\vspace{.2in}
\psadobenomag{task.eps}
\end{minipage}
%\epsfbox{task.eps}
\caption{Compilation of task parallelism}
\label{fig:tasks}
\end{figure}

The compiler can map and schedule instances of
task-subroutines in any way,
as long as  sequential execution results are guaranteed\footnote{Assuming
that input and output parameters are specified correctly.}.
During compilation,
first the input and output directives are analyzed and a task level 
data dependence and communication graph is built.
Next, the task-subroutines are grouped into {\em modules}.
All task-subroutines in the same module are mapped 
to the same set of processors. The modules may be {\em replicated}
to generate multiple module instances, with each instance receiving
and sending data sets in round robin fashion. Constraints on
replication are discussed in section~\ref{sec:replicate}.
Finally, a subset of the
processor array is assigned to each module instance.
All mapping decisions are controlled by the compiler, possibly
using other tools.

Figure~\ref{fig:tasks} shows the application of these steps on
a small example program. The task level dependence graph
is built from the program and then clustered into modules
M1 and M2. Module M2 is replicated to two instances, and the
program is mapped onto the machine.

Our approach to task parallelism is different from
Fortran~M~\cite{FoCh92} in the sense that data and task parallelism
are strongly integrated; both are expressed and compiled uniformly
using a sequential Fortran program and additional directives.  An
effort is underway to integrate task and data parallelism more tightly
in Fortran~M~\cite{chandy93}.  An alternate approach, taken in
Jade~\cite{LaRi91} is to express all parallelism as coarse grain
tasks, and make scheduling decisions at runtime. Most applications
have components that have simple data parallelism, and we believe that
it is extremely important to use an optimizing data parallel compiler
for integrated task and data parallel systems, and statically schedule
and optimize computations and communication whenever feasible.

\section{The mapping problem}
The mapping framework discussed in this paper is applicable to
programs that contain a set of task-subroutines that process a long
stream of inputs in pipelined fashion.  The task-subroutines are
already partitioned by the data parallel compiler to run on a variable
number of processors, and their parallel behavior is characterized by
the following:

\begin{itemize}
\item Response time function: The execution time of
      a task-subroutine call as a function of the number
      of processors.
\item Communication time function: The time to transfer
      data between a pair of task-subroutines as a function
      of the number of processors executing the subroutines,
      and whether or not the pair of task-subroutines is on
      same or different sets of processors.
\item Memory requirement function: The memory requirement
      on each processor for a task-subroutine as a function
      of the number of processors.
\end{itemize}

These characteristics are assumed to be independent of the input data
set, and known at compile time.  They can be determined experimentally
or analytically. In related research, we have shown that a small
number of program executions are sufficient to determine them with
sufficient accuracy for the purpose of identifying a good
mapping~\cite{Subh93b}.

The problem that we are solving in this paper is how to best use the
available processors to maximize the {\em throughput} of the system.
This problem is fundamentally different from the many partitioning
problems addressed in the literature
(e.g.,~\cite{Sark89,Bokh87,nicol91}) for several reasons: (1) the
task-subroutines are to be mapped to a groups of processors, not
individual processors, (2) the computation and communication costs are
functions of the number of processors, not constants, and (3) the
objective is to maximize throughput, not to minimize execution time.

More recently, the problem has been addressed by other researchers,
but not in its general form. Choudhary et. al.~\cite{CNNS94} discuss
the processor allocation problem, but do not address clustering or
replication and do not take into account communication costs that are
not part of the response function.  Cheung and Reeves~\cite{ChRe93}
address replication for a single task-subroutine, but not clustering.
The Paradigm compiler~\cite{RaBa93} attempts to minimize execution
time for a task graph, but the problem is different from maximizing
throughput for a stream of inputs. We also take memory requirements
into account, which are not considered by others, but are of
fundamental importance in using a parallel machine.

It is necessary to address the more general problem to be able to
develop real applications with a parallelizing compiler, even though
the complexity of the problem forces a more heuristic approach.
Computer vision, image processing, and signal processing are examples
of application domains where many programs fit our model, and pure
data parallel computing often leads to disappointing results.  In this
paper we have developed a methodology for mapping task and data
parallel programs that can be used automatically, or by the programmer
with a profiling tool.  We use a set of example programs from these
domains to show how the framework has been used to find good mappings
that combine task and data parallelism.

\section{Characteristics of  parallel programs}
\label{sec:model}
A parallel program in our language contains a set of data
parallel task-subroutines, with communication between them.
We identify three measurable properties of parallel programs
based primarily on communication characteristics
and sizes of data structures.

\begin{enumerate}
\item Scalability
\item Inter-task communication
\item Memory requirement
\end{enumerate}

We now describe  these characteristics, and will use them as the basis for 
our mapping framework discussed in the next section.

\subsection{Scalability}
As the number of processors allocated to a fixed size computation
(or a task) is increased, the average processor
efficiency often decreases. There are a variety
of program characteristics that influence scalability,
particularly the communication overhead and the size
of distributed data structures.
We identify four types of programs and show their
scalability profile in  Figure~\ref{fig:plot}.
Plot A represents a perfectly scalable
parallel program with no communication.
 Plot B represents a  parallel program with
communication, where average processor efficiency decreases
gradually due to increasing  communication overhead,
as the number of processors increases.
Plot C represents a purely sequential 
computation, which cannot use more than one processor. Plot
D represents a computation which scales perfectly up to a
fixed number of processors, but cannot use any more processors.
In general a program represents a mix of these behaviors.

\begin{figure} [p]
\centerline{\psadobenomag{plot.eps}}
\caption{Scalability of computations}
\label{fig:plot}
\end{figure}


\subsection{Communication between task-subroutines}
The volume of data to be transferred in a communication step
between a pair of  task-subroutines is a feature of the program,
but the actual cost of a communication step depends on several
other factors; the number of processors
executing  the task-subroutines, the data distributions,
and the communication style used, which is determined by whether
the two task-subroutines are in the same or different modules.

\subsubsection*{Number of processors}
The cost of
communication depends on the number of processors assigned to
the corresponding task-subroutines. For communicating a fixed amount of data,
the per processor efficiency of communication decreases with
increasing  number of processors, while the total communication time
can increase or decrease.

\subsubsection*{Data distributions}
When the communicating task-subroutines belong to the same module, the
cost varies from zero, when the distribution of data is the same in
both subroutines, to the cost of a fine grained all-to-all
communication, when, for example, one distribution is a transpose of
the other.  In general, the granularity, and hence the overhead of
communication, is dependent on the data distributions.

\subsubsection*{Communication style}
When a pair of task-subroutines is in the same module
(data parallel style), the data exchange is a local redistribution
of data, and
when in different modules (task parallel style), the data is sent
from one set of processors to another set.
These two are different problems, and the cost characteristics can be
different even if the communication patterns are similar.
In the implementation used for our examples,
 communication between modules is done
systolically over a single channel, which
has an extremely low overhead, but limited bandwidth, while  communication
inside a module is done using a message passing library, which has
a relatively higher overhead, but more parallelism and bandwidth.
%The nature of the
%tradeoffs between these styles of communication is illustrated
%in Figure~\ref{fig:plot2}.

%\begin{figure} [p]
%\centerline{\psadobenomag{plot2.eps}}
%\caption{Tradeoff between communication styles}
%\label{fig:plot2}
%\end{figure}

\subsection{Memory requirements}
The memory requirement  of a task-subroutine is an important parameter.
In a multicomputer that has nodes with a fixed amount of uniform
storage, the memory requirement determines the minimum number
of processors needed to execute a task-subroutine, or a module
containing a set of task-subroutines.
In the presence of a memory hierarchy, the memory
requirement  plays an important role in determining overall
performance.

\section{A framework for mapping programs}
\label{sec:framework}
In this section, we outline a framework for reasoning about the
mapping of task and data parallel programs onto multicomputers.
The objective of program mapping is to minimize the execution
time of the whole program, that is to maximize throughput.
 During execution, a processor is in one of the following
states:
\begin{itemize}
\item Executing a task-subroutine 
\item Communicating between task-subroutines
\item Idling due to load imbalance
\end{itemize}

In general, the most important criterion
for high throughput is efficient execution 
of task-subroutines. This may  require that the processors
allocated to individual computations should be small, if
the subroutine does not scale well.
However, reducing the  number of processors assigned to a computation
may not be profitable, or even feasible, due to other factors.
We state the main steps in the process of program mapping, and 
discuss the tradeoffs that have to be made at each step.
 
Following are the major steps in deriving an efficient mapping 
for a parallel program expressed as data parallel task-subroutines:
\begin{description}
\item [Step 1:] Cluster task-subroutines into modules.
\item [Step 2:] Allocate processors to modules.
\item [Step 3:] Replicate modules into module instances.
\end{description}
We present a general framework to make these decisions in 
the above order, although the steps are dependent. In particular,
we have to assess the impact of replication before we make
a clustering decision.

\subsection{Clustering}
The first step is to cluster the task graph obtained from
the program into modules. The nodes in the task
graph are data parallel procedures which can be potentially
replicated, and whose cost is a function of the number of 
processors assigned to them. The edges correspond to
communication cost, which varies with the number of processors,
and also depends on whether the the corresponding nodes are assigned
to the same or different modules.

We first assign a module to each task-subroutine, and compare pairs
of neighboring modules to determine if they should be merged into a
single module. The order in which pairs of modules are tested for
potential merger is determined heuristically based on the following
priority scheme:
\begin{itemize}
\item Less scalable modules are prioritized, since their 
mapping is most important for efficient execution of the program.
\item A pair of modules that share most of the data structures
are prioritized, since a decision to combine them has minimal
impact on other decisions.
\end{itemize}
We now discuss the procedure for deciding
if a pair of modules should be merged or not.
We first examine the considerations for this decision:

\begin{description}
\item[Replicability:] If one of the modules cannot be
replicated, merging them implies that neither of them can be 
replicated, hence the unmerged modules only will benefit
from replication.
\item[Memory constraints:] Merging the modules increases the
memory requirement, and hence the minimum number of processors
that the merged module needs to execute can be higher than
the minimum  number of processors needed by separate modules.
This can restrict (or prevent) replication,
and  limit the best  performance, as
illustrated in Figure~\ref{fig:plot}.
\item[Inter-task communication cost:] If the task-subroutines
are merged, the communication step is a data redistribution,
and can be zero if the data structures have the same 
distribution in both task-subroutines. If they are not merged,
the data has to be transferred between sets of
processors.
\end{description}

We now state the specific criterion for making this decision
assuming replication is possible for both modules. The obvious 
approach is to estimate the performance in the two situations, but
the difficulty is that the performance of non-scalable modules 
depends on the number of
processors assigned to each module, and that is not known. 
However,  when the modules can be replicated, an optimistic
estimate of performance can be obtained by assuming that  
each instance of the module will be assigned the minimum number 
of processors that are necessary due to memory or I/O requirements.
So, if we are considering merging modules $M_1$ and $M_2$ into a single
module $M_{12}$, and the {\em minimum} number of processors required to
fit the data sets of corresponding modules is $P_1$, $P_2$ and $P_{12}$,
we use these processor numbers to estimate computation and communication 
performance.

\subsection{Processor allocation}
The allocation of processor is based on the execution time estimates
of the modules discussed in the last section.
Methods discussed in~\cite{CNNS94} can be used when fixed
estimates are not possible.

\subsection{Replication} \label{sec:replicate}
Replication improves the performance by reducing the
execution overhead of modules, as well as
the cost of communication  between them, since each
instance executes on a smaller number of processors.
The effect is most pronounced for modules that do not
scale well. However, it is not possible to arbitrarily
replicate modules into smaller module instances.
The legality and the degree of maximum possible replication
is determined by the following constraints:
\begin{description}
\item[Dependence] Any module that has carried dependence
(or carries {\em state}) cannot be replicated.
\item[Memory Requirement]
Every module needs a minimum number of processors to
hold the dataset, which bounds the number of module
instances possible for a fixed set of processors.
\item[Latency]
The latency between receiving an input and  generating an output
can increase with replication, even if the throughput increases.
Although we have assumed that our objective is to maximize
throughput, the extent of latency that can be tolerated may be  limited in
real-time applications.
\item[Resource constraints]
Supporting a finer granularity of task parallelism obtained by
replication increases the demand for machine resources. This
can constrain replication either because resources are not
available, (e.g. communication channels) or increase the 
overhead due to congestion. Communication resources
also restrict the number of modules that can interface with
external I/O agents.
\end{description}

For the purpose of determining a clustering, we assume maximum
replication  allowed by the constraints stated above, but
do not consider the effect of resource consumption of one
module on the availability of resources for others. This is
done to obtain  computation and communication time
estimates before clustering into modules is complete.

The actual replication is done in the order of increasing
scalability, that is, the least scalable modules are replicated
first. For each module, we pick the minimum number of processors
that are needed to obtain nearly linear speedup, and replicate
until that number of processors per instance is reached, or
one of the constraints prevents further replication. The process
stops when all the modules are processed, or one of the essential
resources is exhausted.

\section{Example programs}
In this section, we apply the framework developed in the previous
section to three realistic task and data parallel programs: fast
Fourier transform, narrowband tracking
radar~\cite{SGMR92}, and multibaseline stereo imaging~\cite{webb93b}.
We briefly describe each example program, show how it is implemented as an
Fx program, and show how to use the framework to reason about mapping
task parallel programs. We present performance results on a 64
processor iWarp array~\cite{BCCG90} to validate our approach.  

Each example program has unique properties that exercise different
parts of the mapping framework.  Each program operates on a stream of
input, and generates a stream of output.  While the I/O issues
involved are beyond the scope of this paper, we assume that the I/O
streams can be partitioned into up to four interleaved I/O substreams.

\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 $y = F_n
x$ where $x$ and $y$ are complex vectors, and $F_n = (f_{pq})$ is an
$n \times n$ matrix such that $f_{pq} = \omega_n^{pq}$ where $\omega_n
= \cos(2\pi/n) - i \sin(2\pi/n) = e ^ {-2 \pi i /n }$ and $i =
\sqrt{-1}$.  An FFT algorithm exploits structure in $F_n$ to compute
the DFT using $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.

Large 1D~FFTs, 2D~FFTs, as well as higher dimensional FFTs can be
computed as a sequence of small, independent, data parallel  
1D~FFTs with transpose-like operations between them. 
If $n=n_1 n_2$, then a 1D~FFT can be computed using a collection of
smaller independent 1D~FFTs~\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$. The final result vector $y$
is obtained 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$, we can compute a 2D~FFT by
performing $n_1$ $n_2$-point 1D~FFTs on the rows of A, followed by
$n_2$ $n_1$-point FFTs on the columns of $A$.

The issues involved in mapping all FFTs are very similar.
Here we choose a $128 \times 128$ 2D~FFT as an example.
An Fx 2D~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 rffts(A)\\
C\$\>\>     {\bf output:} A\\
\>\>       call cffts(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.  Figure
\ref{fig:ffttask} shows the task graph for one iteration of the
parallel section.  Subroutine {\tt rffts} 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 cffts} 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.

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

\subsection*{Using the framework to map the FFT program}
The 2D~FFT program consists of two data parallel task-subroutines
({\tt rffts} and {\tt cffts}) separated by a transpose-like inter-task
communication step that copies an array distributed by rows to an
array distributed by columns.  Here then the steps we take to map
a $128 \times 128$ 2D~FFT onto a 64-processor system, according to the
framework described in Section~\ref{sec:framework}:

{\em Step 1: Clustering.} There are only two possible clusterings:
assigning the two task-subroutines to the same module or to different
modules. For the FFT program, both the task-subroutines can be 
replicated  because the input data sets are independent
and because neither of the task-subroutines carries any state.
The maximum degree of replication for either task-subroutine
is 4, since there are a maximum of four I/O streams (as discussed in the
beginning of the section).
According to section~\ref{sec:framework}, we have to compare
performance assuming maximum replication, to make the clustering decision.

If the two task-subroutines are in the same module, then after maximal
replication, each instance of the module
executes on 16 processors, and the inter-task communication step uses
message passing to redistribute a $128 \times 128$ complex array on 16
processors. From  performance estimates,
one iteration of {\tt cffts} and {\tt rffts} executes  in
169.7 ms and 166.0 ms respectively, and the communication step takes 81.0 ms, 
yielding a total execution time of 416.7 ms. for one iteration
on 16 processors. 

On the other hand, if the two task-subroutines are in
different modules, then after maximal replication, each instance of a
module executes on 8 processors, and the communication step is a
transfer of a $128 \times 128$ complex array between two sets
of 8 processors. From the performance estimates,
one iteration of {\tt cffts} and {\tt rffts} executes  in
340.4 ms and 324.1 ms respectively on 8 processors,
 and the communication step takes 55.1 ms.
The average execution time for one iteration on per 16 processors
(taking pipelining into account and assuming perfect load balancing)
is ((340.4 + 55.1)*8 + (324.1+55.1)*8)/16 = 387.3 ms.

Since keeping the task-subroutines separately gives a slightly better
relative performance estimate (387.3 ms versus 416.7 ms), we 
put them in different modules. Note that the modules are computationally
scalable, and the decision was mainly based on communication estimates.

{\em Step 2: Processor Allocation.} Since both modules perform nearly
identical computations, each is
allocated half (32) of the total number (64) of processors.
  
{\em Step 3: Replication.} Replication is desirable
for the FFT because it reduces the the total execution time by
reducing the cost of the inter-task
communication step, but is limited to 4 due to I/O constraints.
So we replicate each module 4 times, with each instance of a module
running on 8 processors. 
The resulting mapping for the $128 \times 128$ 2D~FFT program is shown
in Figure~\ref{fig:ffttask}.  

Figure \ref{fig:fftnumbers1} shows the speedups of the
different mappings of the FFT program
relative to a simple data parallel mapping, for different problem sizes.
For the $128 \times 128$ FFT, the mapping we derived with the
framework outperforms all others(albeit very slightly over the
``1 module replicated 4 ways'' mapping). Similarly, we were able to
derive the best mappings for FFTs of other sizes, but the discussion
is omitted for brevity. Notice that the optimal clustering depends
on the problem size, but a higher degree of replication
always improves performance.

\begin{figure}[p]
\begin{center}
\begin{tabular}{|l|c|c|c|c|c|c|}
\hline
&\multicolumn{3}{c|}{one module} &  \multicolumn{3}{c|}{two modules} \\ 
\cline{2-7} 
size     & pure data   & 2-way  & 4-way & no     & 2-way       & 4-way \\
 & parallel & replication & replication  & replication & replication & replication\\
\hline
$64 \times 64$   & 1 & 1.88 & 2.61 & 2.08 & 2.65 & 3.04 \\
$128 \times 128$ & 1 & 1.15 & 1.28 & .91 & 1.12 & 1.29\\
$256 \times 256$ & 1 & 1.09 & 1.16 & .79 & .94 & 1.08 \\
\hline
\end{tabular}
\end{center}
\caption{Measured speedup of FFT relative to a pure data parallel mapping.}
\label{fig:fftnumbers1}
\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{SGMR92}. It is a
particularly interesting benchmark for studying task parallelism
because of its hard real-time requirements, and because the size of
the input data set is limited by physical properties of the radar
sensor.  The amount of available low-level data parallelism is
limited, so additional parallelism must come from higher-level task
parallelism.

The radar 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.

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 dgen(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
Like FFT, the program consists of a parallel section with calls to two
task-subroutines inside a loop that iterates $m$ times.  Figure
\ref{fig:radartask} shows the task graph. Task-subroutine {\tt
dgen} acquires the data from all 4 channels and sends it to
task-subroutine {\tt compute}, a data parallel routine that performs
steps (1)--(4) above.  The data parallelism is in the form of a
parallel loop where each loop iteration operates on a single column of
the corner-turned data set.

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

\subsection*{Using the framework to map the radar program}
The radar program has two characteristics that dominate the choice of
a mapping. First, the data generation task-subroutine must run on a
single processor because it acquires the input data from a sensor that
is attached to a particular processor. Second, since there are only 10
columns in the input array after it has been corner turned, and an FFT
computation is done on each column, the amount of data parallelism
available to the {\tt compute} task-subroutine is quite small.
Following are the steps we take to map the radar program onto 64
processors.

{\em Step 1: Clustering.} As with the FFT example, there are two
task-subroutines and two possible clusterings.
Replication of the {\tt dgen} task-subroutine is not
possible because it must run on one processor. However, replication of
the {\tt compute} task-subroutine is both possible and desirable;
possible because the input data sets are independent, desirable
because each instance of the {\em compute} task-subroutine is able to
efficiently use only a small number of processors.  If we assign the
task-subroutines to the same module, then we eliminate the possibility
of replicating the {\tt compute} task-subroutine.  So we assign the
task-subroutines to different modules, as shown in
Figure~\ref{fig:radartask}.

{\em Step 2: Processor Allocation.} For the radar program, this step
is trivial.  We assign one processor to the {\tt dgen} module and the
remaining 63 processors to the compute module.

{\em Step 3: Replication.} Since the input data set is small,
and {\tt compute} is not scalable beyond 10 processors, replication
is desirable. However,
latency constraints in the specification of
the application as well as resource constraints in the communication
system of the target system limit us to 4 replicated modules, each
running on 10 processors, the maximum number that a module can use 
efficiently in this implementation.

The resulting mapping is shown in Figure~\ref{fig:radartask}.
Figure~\ref{fig:radarnumbers} gives the measured performance of the Fx
radar program when compiled with different degrees of replication.
The near linear speedups illustrate the value of replication, for
programs like the radar program that operate on small data sets.
\begin{figure}
\begin{center}
\begin{tabular}{|c|c|c|c|}
\hline
replication & 1-way & 2-way & 4-way\\
\hline
\hline
time/input(ms) & 49.6 & 24.8 & 12.4 \\
\hline
\end{tabular}
\end{center}
\caption{Performance of radar under different degrees of replication.}
\label{fig:radarnumbers}
\end{figure}

\subsection{Multibaseline stereo}
The multibaseline stereo uses an algorithm developed at Carnegie
Mellon that gives greater accuracy in depth through the use of more
than two cameras~\cite{okutomi93}. It is an interesting program for
studying task parallelism because it contains significant amounts of
both inter-task and intra-task communication\cite{webb93b}, and
the size of data sets is fixed.
Our implementation is adapted from a
previous data-parallel implementation written in a specialized image
processing language~\cite{webb93a}.

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 = 0,\ldots,15$, 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 corresponding 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.  Finally, the depth of each pixel is
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 dgen(R,M1,M2)\\
C\$\>\>    {\bf output:} R,M1,M2\\
\>\>       do d = 0,15\\
\>\>\>		call diff(R,M1,M2,DIFF,d)\\
C\$\>\>\>       {\bf input:} R,M1,M2\\
C\$\>\>\>       {\bf output:} DIFF\\
\>\>\>		call error(DIFF,ERR(:,:,d),d)\\
C\$\>\>\>       {\bf input:} DIFF\\
C\$\>\>\>       {\bf output:} ERR(:,:,d)\\
\>\>       enddo\\
\>\>       call min(ERR,DISP)\\
C\$\>\>    {\bf input:} ERR\\
C\$\>\>    {\bf output:} DISP\\
\>     enddo\\
C\$\>   {\bf end parallel}
\end{tabbing}
\rm
\normalsize
Figure \ref{fig:stereotask} shows the task graph.  Task-subroutine
{\tt dgen} acquires three $256 \times 240$ images from the cameras.
Each of the 16 instances of the {\tt diff} task-subroutine is a
perfectly data-parallel routine that converts the three 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 its neighbors before the error
image can be produced. The outputs from the {\tt error} tasks are
passed to {\tt min}, which applies a {\em min}-reduction to produce
the disparity image, and then displays the corresponding depth image.

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

\subsection*{Using the framework to map the stereo program}
The relative performance of different mappings of the stereo program
is determined primarily by (1) the cost of the inter-task
communication steps that pass images from one task-subroutine to the
next (without changing the distribution), and (2) the cost of the
intra-task communication steps dominated by the {\tt error} 
task-subroutines.  Following are the steps we take to
map the stereo program onto 64 processors:

{\em Step 1: Clustering.}
For all pairs of
communicating task-subroutines in the stereo program, there is no
inter-task communication cost if they are mapped to the same module,
since the data distributions are the same.  However, if we place any
two communicating task-subroutines in different modules, there is a
significant data movement and communication cost, and there is no gain
in terms of replicability or memory requirements.  So we decide to
place all task-subroutines in a single module, as shown in
Figure~\ref{fig:stereotask}, and use replication to enhance
performance.

{\em Step 2: Processor Allocation.} The clustering from Step 1 consists
of one module, to which we assign all 64 processors.

{\em Step 3: Replication.} Memory constraints in our target system
limit the degree of replication to four modules, each running on 16
processors. The final mapping is shown in Figure~\ref{fig:stereotask}.
Figure~\ref{fig:stereonumbers1} shows the measured performance of the
Fx stereo program compiled as 1 (i.e. purely data parallel),
 2, and 4 replicated modules. The
higher throughput of the 4-way replicated case validates the decision
to use replication.  However, notice that
while a 4-way replication roughly doubles the throughput, it roughly
doubles the latency too. Depending on the requirements of a particular
application of the stereo program, this may or may not be a reasonable
tradeoff. A system striving to minimize latency would potentially
arrive at a different mapping.

\begin{figure}[p]
\begin{center}
\begin{tabular}{|c|c|c|}
\hline
degree of   & latency           & throughput\\
replication & msec/frame & frames/sec\\
\hline
\hline
1 & 112 & 8.9 \\
2 & 155 & 12.9 \\
4 & 231 & 17.2 \\
\hline
\end{tabular}
\end{center}
\caption{Performance of $256 \times 240$ stereo with replication }
\label{fig:stereonumbers1}
\end{figure}

\comment{
\subsubsection*{Latency and throughput tradeoffs in the stereo program}
The stereo
program captures some potentially important tradeoffs between latency
and throughput.  Latency is the
time required to produce a given frame (an output depth image) from
the time the corresponding trio of images is received from the camera.
Throughput is the aggregate number of frames produced each second.
The mapping procedure described above
attempted to optimize the throughput, resulting
in a 4 way replicated mapping.
Figure~\ref{fig:stereonumbers1} shows the measured performance of
the Fx stereo program compiled as 1, 2, and 4 replicated
modules. Notice that a 4 way replication roughly doubles the
throughput, but roughly doubles the latency too. Depending
on the requirements of a particular application of the stereo program,
this may or may not be a reasonable tradeoff. For example, 
the mapping shown in Figure~\ref{fig:stereotask2} can potentially 
improve the throughput and latency in the program by speeding up
processing of each data set, rather than processing multiple 
data sets simultaneously. Further discussion on improving latency
is beyond the scope of this paper and we merely wish to point out
that optimization for latency and throughput are very different
problems, and one or the other is relevant depending on the 
situation.

\begin{figure}
%\epsfxsize=4in
\centerline{\epsfbox{stereotask2.eps}}
\caption{Stereo task graph clustered into 4 modules}
\label{fig:stereotask2}
\end{figure}
}
 
\subsection{Summary of example programs}
For each example program, we are able to use the framework to
derive mappings that exploit a mix of task and data parallelism to
achieve better performance than that obtained from a straightforward
data parallel mapping.  Each example program exercises different
aspects of the framework. In the FFT program, the parallelization
overhead is dominated by the communication step between the two FFTs.
By clustering into two modules and using replication,  the efficiency
of communication is maximized.  In the narrowband tracking radar program, the
main computation routine can effectively use only 10 processors, and
improvement is possible by replicating this routine, which in
turn is made possible by placing the {\em datagen} and {\em
compute} routines into separate modules.  In the multibaseline stereo
program, the cost of communication inside task-subroutines is reduced
by replicating the entire computation,  which results in each
task-subroutine instance executing on fewer processors with better
locality. 

\section{Automatic mapping with real world constraints}
In this paper we have demonstrated that a variety of applications can
benefit from using a mix of task and data parallelism, and that good
mapping decisions can be based on well understood program properties.
In related work, we have used this analysis to develop an automatic
tool that obtains profile data from program executions, and uses it to
automatically generate a good mapping~\cite{Subh93b}.

Our experience is that real applications in signal and image 
processing often have constraints that cannot be captured easily by
a fully automatic system, hence user involvement in the mapping
process is important. There are several reasons for this. The
mapping is affected by the supporting I/O framework which is 
hard to model in general. For instance, in the stereo program
only specific processors are connected to the camera input, and
the number of I/O streams is determined  externally.
In the radar program, only a single processor receives the input
stream, and hence the corresponding module cannot be replicated.
Another factor that is difficult to model automatically is the
memory hierarchy. While we are able to measure the total memory 
requirements of a program, it is difficult to model the cache
and secondary memory behavior. Also, the effects of congestion
in the communication system are extremely difficult to model in
a general framework. Finally, some scalability characteristics
can be difficult to derive automatically from measurements.
For example in the radar program, the program performance 
scales linearly up to 10 processors, but not at all after that,
due to the parallel dimension of the main structure being of
size 10.

In summary, we believe that various automatic tools can be invaluable
to a programmer for deciding a mapping, but the current state of the
art does not allow a satisfactory fully automatic approach.  With the
Fx toolchain, the programmers can use an automatic tool to guide
mapping, or derive a mapping using a profiling tool and the framework
presented in this paper. In this paper we have focused on the basic
methodology for obtaining a good mapping which can be used in
different ways.


\section{Concluding remarks}
The performance of any parallel program depends on how it is mapped onto a
parallel machine.  For many programs, a simple strategy that focuses only
on task or data parallelism is not adequate, and better results can be
obtained using both task and data parallelism.

Finding good mappings of parallel programs is one of the key problems
in exploiting task and data parallelism together.  This paper
identifies the fundamental program properties that define feasible
mappings and determine their performance. The communication
requirements of a parallel program determine the relative performance
of different mappings. The memory requirements determine if a mapping
is feasible, since the local memory of each node is finite.  Based on
the analysis of these properties and machine parameters, we have
developed a framework for finding good mappings.

We have demonstrated the mapping process using three realistic
parallel programs developed with our compiler. In all cases, the
mapping obtained was more efficient than a simple data parallel
mapping.

\bibliography{/afs/cs/project/iwarp/member/jass/bib/all,/afs/cs/project/iwarp/member/droh/bib/defs,/afs/cs/project/iwarp/member/droh/bib/me,/afs/cs/project/iwarp/member/droh/bib/refs}

%/afs/cs/project/iwarp/member/droh/bib/compiler
\end{document}

