% Technical report on Global Arrays and High-Performance Fortran
% prepared for Battelle/Pacific Northwest National Laboratory under
% contract 313130-A9E
%
% $Id: report.tex,v 1.1 1999/04/01 22:02:35 dbc Exp $
%

\documentclass[12pt]{article}

\usepackage{psfig}

\begin{document}

\title{Interfacing Global Arrays and ARMCI with the PCRC
library, Adlib\thanks{Prepared for Battelle/Pacific
Northwest National Laboratory under contract 313130-A9E.}}

\author{Bryan Carpenter \\
Northeast Parallel Architectures Center\\
Syracuse University\\
Syracuse, NY 13422-410 \\
{\em dbc@npac.syr.edu} \\
\\
Jarek Nieplocha \\
Pacific Northwest National Laboratory, \\
Richland, WA 99352 \\
{\em j\_nieplocha@pnl.gov}
}


\maketitle

\begin{abstract}
This document reports work undertaken at NPAC, Syracuse under
the DOE Global Array Extension Project.  This work was intended to
investigate the feasibility of interfacing, and perhaps eventually
integrating, GA, ARMCI and the Parallel Compiler Runtime
Consortium library, Adlib.  In particular, we have
reimplented parts of the Adlib library in terms of ARMCI, and also
produced a version of GA which internally uses an Adlib-compatible array
descriptor.
\end{abstract}

\tableofcontents

\newpage

\section{Introduction}

\subsection{Context}

The Department of Energy Advanced Computational Testing and Simulation
(ACTS) project combines multiple programming tools and environments to
support application codes addressing complex multidisciplinary
problems.  The Global Array (GA) library developed at the Pacific
Northwest National Laboratory (PNNL) is one of the ACTS components.  In
ACTS the GA tool kit will be extended to provide support for higher
dimensional arrays and additional distribution formats.
The Adlib \cite{NPAC_PCRC_kernel} parallel run-time
library has capabilities that can be used to implement these new features.

Development of the Adlib library was completed at NPAC, Syracuse in the
framework of the three-year Parallel Compiler Runtime Consortium (PCRC)
project.  It incorporates experience gained over the years in parallel
distributed array technology.  PNNL determined that it would be
cost effective to transfer some of the multidimensional array
capabilities from Adlib to Global Arrays, rather than to develop them
from scratch at PNNL.

PNNL proposed to build an implementation of the GA tool kit
Application Programming Interface (API) on top of a version of the PCRC
library, augmented with the ARMCI aggregate remote memory copy
library.  ARMCI is being developed by PNNL.  Towards this goal an
initial study phase in fiscal year 1998 (reported on here) aimed to
establish the basic feasibilty of this approach.  We had to provide
evidence that a new layered implementation would not lead to an
unacceptable degradation of performance.  For example, it had to be
demonstrated that the global-to-local subscript conversion operations
supported by the Adlib array descriptor (suitably optimized if
necessary) would not compromise the low latency of the current {\tt
ga\_get}/{\tt put} operations.  It was also considered that it would be
very useful to reimplement part of the Adlib collective communication
library on top of ARMCI, and compare the remote-memory-access approach
with the current message passing implementation.

\subsection{Summary of achievements}

In the pilot study parts of the Adlib library were reimplemented on top
of the ARMCI library.  Early benchmarks show that the new
implementation is significantly faster than the original MPI
implementation on shared memory platforms.  In a second part of the
study, a version of Global Arrays was created that internally uses the
Adlib distributed array descriptor to parametrize the decomposition of
multidimensional global arrays, and uses associated functions from
Adlib to perform address translations.  The modified GA was
benchmarked.  The results indicate that the overhead of introducing the
object-oriented DAD of Adlib into the GA implementation are
small---probably less than a microsecond for a {\tt ga\_get}
operation.


\subsection{This report}

Section \ref{sec:background} reviews relevant aspects of the Adlib
runtime library.  For completeness it also briefly reviews the
Global Arrays toolkit and ARMCI.  Section \ref{sec:AdlibOnARMCI}
describes work on using ARMCI in the implementation of Adlib.
In particular it describes how various existing Adlib communication
schedules were reimplemented in terms of ARMCI, and gives benchmark
results.  It also discusses how ARMCI will allow the Adlib API
to be extended with functions for one-sided-communication. 
Section \ref{sec:GAonDAD} describes the pilot work on using
Adlib-based techniques in the implementation of GA.

\section{Background\label{sec:background}}

\subsection{The Adlib library}

The Adlib library was completed in the Parallel Compiler Runtime
Consortium project \cite{Common_runtime} project.  It is a high-level
runtime library designed to support translation of data-parallel
languages \cite{NPAC_PCRC_kernel}.  Initial emphasis was on High
Performance Fortran (HPF), and two experimental HPF translators used the
library to manage their communications \cite{Zhang97,merlin96}.
Currently the library is being used in the {\em HPspmd} project at NPAC
\cite{Carpenter98}.  It incorporates a built-in
representation of a distributed array, and a library of communication
and arithmetic operations acting on these arrays.  The array model is
more general than GA, supporting general HPF-like distribution formats,
and arbitrary regular sections.  The existing Adlib communication
library emphasizes collective communication rather than one-sided
communication.

\subsubsection*{Array descriptor\label{sec:DAD}}

An important part of the Adlib library is its characteristic {\em
Distributed Array Descriptor} (DAD).  Adlib is implemented in C++, and
the underlying array descriptor is implemented as an object of type {\tt
DAD}.  A DAD object describes how elements of a particular array are
distributed across available processors.  It has a simple structure
(Figure \ref{fig:DAD})
containing three fields: an integer value ($r$, say) defining the {\em rank}
(dimensionality) of the array, a process {\em group} object, and a
vector of $r$ {\em map} objects, one for each dimension of the array.
The group defines a
multidimensional process grid embedded in the set of processes
executing the program, or some subset.
Each map object consists of an integer local {\em memory stride} and
a {\em range} object.

\begin{figure}[tb]
  \begin{center}
    \leavevmode
    \psfig{figure=DAD.eps,height=2.6in}
    \caption{Structure of the Adlib DAD for a 3-dimensional array}
    \label{fig:DAD}
  \end{center}
\end{figure}

A range object describes the {\em extent} (size) and distribution
format of one dimension of a distributed array.  Adlib defines a class
hierarchy of different kinds of range object (Figure \ref{fig:ranges}).
Each subclass represents a different kind of distribution format for an
array dimension.  The simplest distribution format is {\em collapsed}
(sequential) format in which the whole of the array dimension is mapped
to the local process.  Other distribution formats (motivated by High
Performance Fortran) include {\em regular block} decomposition, {\em
simple cyclic} decomposition, and {\em block cyclic} decomposition.  In
these cases the index range (thus array dimension) is distributed over one of
the dimensions of the process grid defined by the group object in the
DAD.  All ranges in a specific DAD must be distributed over different
dimensions of this grid, and if a particular dimension of the grid is
targetted by none of the ranges, the array is said to be {\em
replicated} in that dimension\footnote{So there is no direct relation
between the array rank and the dimension of the process grid:
collapsed ranges means the array rank can be higher; replication allows
it to be lower.}.  Some of the range classes allow
{\em ghost extensions} to support stencil-based computations.

\begin{figure}[tb]
  \begin{center}
    \leavevmode
    \psfig{figure=ranges.eps,height=2.2in}
    \caption{The Adlib {\tt Range} hierarchy}
    \label{fig:ranges}
  \end{center}
\end{figure}

The factorization of the DAD into separate map and range objects for each
dimension is convenient for several reasons.  For example it simplifies
construction of Fortran-style {\em regular sections} of arrays.  It
means that arrays with common distribution patterns for selected
dimensions can share the same range objects for those dimensions.
Also, the range object is a natural locus for the methods used to
enumerate locally held index blocks for a distributed array dimension.
These methods are used in translation of {\em forall} loops and similar
data parallel constructs, and inside collective communication
libraries.

The DAD itself does not contain the address of the memory segment where
locally held elements are stored.  Nor does it contain any
information about the type of those elements.  This information must be
maintained by the user of the kernel library or (more often) in
higher-level interface code wrapped around the kernel library calls.  This
arrangement allows the kernel DAD to be independent of details of
the user-level programming language and the memory managment scheme
used to allocate local array elements.
As an example of the kind of interface that is possible, there
is a C++ user-level interface to Adlib called
{\em ad++}.  This defines a series of template classes
\small
\begin{verbatim}
        template<class T>
          class Array1 ;

        template<class T>
          class Array2 ;

        template<class T>
          class Array3 ;

        ...
\end{verbatim}
\normalsize
These implement distributed container classes for elements of type {\tt
T} (arrays of rank 1, 2, 3, \ldots).  They are derived from the base class
{\tt DAD}, but add a field pointing to vectors of local elements; their
constructors allocate these vectors.

\subsubsection*{Communication schedules\label{sec:schedules}}

The Adlib communication library supports {\em collective operations}
on distrib\-uted arrays.  A call to a collective operation must be
invoked simultaneously by all members of some {\em active process
group}, which may or may not be the entire set of processes executing
the program.

Communication patterns supported include HPF/F90 intrinsics such as
CSHIFT and TRANSPOSE.  More importantly they include the
regular-section copy operation, {\tt remap}, which copies elements
between shape-con\-forming array sections regardless of source and
destination mapping.  Another function, {\tt writeHalo}, updates ghost
areas of a distributed array.  Various collective gather and scatter
operations allow irregular patterns of access.  The library also
provides essentially all F90 arithmetic transformational functions on
distributed arrays and various additional HPF library functions.

All collective operations in the library are based on communication
{\em schedule} objects.  Each kind of operation has an associated class
of schedules.  Particular instances of these schedules, involving
particular data arrays and other parameters, are created by the class
constructors.  Executing a schedule initiates the communications required
to effect the operation.  A single schedule may be executed many
times, repeating the same communication pattern.  In this way,
especially for iterative programs, the cost of computations and
negotiations involved in constructing a schedule can often be amortized
over many executions.  This paradigm was pioneered in the CHAOS/PARTI
libraries \cite{CHAOS}.  If a communication pattern is to be executed
only once, simple wrapper functions can be made available to construct
a schedule, execute it, then destroy it.  The overhead of creating the
schedule is essentially unavoidable, because even in the single-use
case individual data movements generally have to be sorted and
aggregated, for efficiency.  The associated data structures are just
those associated with schedule construction.

A characteristic example from the standard Adlib communication library is
the remap schedule.  The associated class has public interface:
\small
\begin{verbatim}
  class Remap {
  public :
    Remap(const DAD* dst, const DAD* src, const int len) ;

    void execute(void* dstDat, void* srcDat) ;
  } ;
\end{verbatim}
\normalsize
The constructor is passed the DADs for source and destination
arrays (or sections) and the size of the array element in bytes.
The {\tt execute} member is passed the base addresses for the locally
held elements (the destructor for the schedule is not
displayed here).  The ad++
wrapper function for single use of this schedule with two dimensional
arrays is
\small
\begin{verbatim}
 template<class T>
 inline void remap(const Section2<T>& dst, const Section2<T>& src) {
   Remap schedule(&dst, &src, sizeof(T)) ;

   schedule.execute(dst.dat(), src.dat()) ;
 }
\end{verbatim}
\normalsize
The inquiry {\tt dat} is a member of the ad++ array classes which
returns a pointer to the local data vector.  The {\em remap}
function copies all elements of the source array (or section) to
the destination array (or section)\footnote{The {\tt Section}$R${\tt<T>}
classes are base classes of the {\tt Array}$R${\tt<T>} classes.}.

Top-level schedules such as {\tt Remap}, which deal explicitly with
distributed arrays, are implemented in terms of some
lower-level schedules that simply operate on blocks
and words of data.  These lower-level schedules do not engage in
operations on the DAD or its range and group components.
The lower level schedules are tabulated in Table \ref{tab:schedules}.
Here ``words'' are contiguous memory blocks of constant (for
a given schedule instance) size.  ``Blocks'' are multidimensional
($r$-dimensional) local array sections, parametrized by a vector
or $r$ extents and a vector of $r$ memory strides.
The point-to-point schedules are used to implement collective
operations that are deterministic in the sense that both sender
and receiver have a priori knowledge of all required communications. 
Hence {\tt Remap} and other regular communications such as {\tt Shift}
are implemented on top of {\tt BlockMessSchedule}.  The ``remote access''
schedules are used to implement operations where one side must
inform the other end that a communication is needed.  These negotiations
occur at schedule-construction time.  Irregular communication
operations such as collective {\tt Gather} and {\tt Scatter} are
implemented on these schedules.  The tree schedules are used
for various sorts of broadcast, multicast, synchronization, and
reduction.

\begin{table}[tb]
\caption{Low-level Adlib schedules\label{tab:schedules}}
\vspace{0.1in}
\centerline{
\begin{tabular}{| c | c | c |}
\hline
                & operations on ``words'' & operations on ``blocks'' \\
\hline
Point-to-point  & {\tt MessSchedule }  & {\tt BlockMessSchedule } \\
\hline
Remote access   & {\tt DataSchedule }  & {\tt BlockDataSchedule } \\
\hline
                & {\tt TreeSchedule }  & {\tt BlockTreeSchedule } \\
Tree operations & {\tt RedxSchedule }  & {\tt BlockRedxSchedule } \\
                & {\tt Redx2Schedule } & {\tt BlockRedx2Schedule } \\
\hline
\end{tabular}
}
\end{table}

The original Adlib is implemented on top of MPI.  All MPI
communications are isolated in the low-level schedules.  Ideally, one
might hope that retargetting suitable parts of Adlib to operate on top of
ARMCI would be simply a matter of reimplementing the low-level
schedules.  In practise, as explained in section
\ref{sec:AdlibOnARMCI}, the layered implementation eased the task of
retargeting to ARMCI, but the API of the low-level schedules had to
be changed.

\subsection{ARMCI and Global Arrays}

ARMCI (Aggregate Remote Memory Copy Interface) is a new portable remote
memory copy library developed at PNNL for optimized communication
in the context of distributed arrays.  ARMCI aims to be fully portable
and compatible with message-passing libraries such as MPI or PVM.
Unlike most existing similar facilities, such as Cray SHMEM \cite{Barriuso94}
or IBM LAPI \cite{LAPI}, it focuses on the noncontiguous data
transfers.  ARMCI offers a simpler and lower-level model of one-sided
communication than MPI-2 \cite{MPI2} (no epochs, windows, datatypes,
Fortran-77 API, complicated progress rules, etc.) and targets
a different user audience.  In particular, ARMCI is meant to be used
primarily by library implementors rather than application developers.
Examples of libraries that ARMCI is aimed at include Global Arrays,
P++/Overture \cite{Quinlan94,Bassetti98}, and Adlib.

ARMCI supports three classes of operation
\begin{itemize}
\item
data transfer operations including put, get and accumulate.
\item
synchronization operations---local and global fence and atomic
read-modify-write, and
\item
utility operations for allocation and deallocation of memory
(as a convenience to the user) and error handling.
\end{itemize}
The data transfer operations are available with two noncontiguous
data format:
\begin{enumerate}
\item
{\em Generalized I/O vector.}  This is the most general format intended
for multiple sets of equally-sized data segments, moved between arbitrary
local and remote memory locations.  It extends the format used in the
UNIX {\em readv/writev} operations by minimizing storage requirements
in cases when multiple data segments have the same size.
The associated functions are {\tt ARMCI\_PutV} and {\tt ARMCI\_GetV}.
\item
{\em Strided.}  This format is an optimization of the generalized I/O
vector format for sections of dense multidimensional arrays.  Instead
of including addresses for all the segments, it specifies for source
and destination only the address of the first segment in the set.  The
addresses of the other segments are computing using the stride
information.  The associated functions, {\tt ARMCI\_PutS} and
{\tt ARMCI\_GetS}, take short {\em vectors} of strides as arguments
allowing {\em multidimensional} array segments to be handled in
a single operation.
\end{enumerate}

The Global Arrays (GA) \cite{GA} toolkit provides a shared-memory
programming model in the context of 2-dimensional distributed arrays.
GA has been the primary programming model for numerous applications, some
as big as 500,000 lines of code, in quantum chemistry, molecular dynamics,
financial calculations and other areas.  The toolkit was originally
implemented directly on top of system-specific communication mechanisms
(NX {\em hrecv}, MPL {\em rcvncall}, SHMEM, SGI {\em arena}, etc.).
The original one-sided communication engine of GA had been closely
tailored to the two-dimensional array semantics supported by the toolkit.
This specificity hampered extensibility of the GA implementation.  It
became clear that a separation of the communication layer from the
distributed array infrastructure is a much better approach.  This was
accomplished by restructuring the GA toolkit to use ARMCI hence making
the GA implementation itself fully platform independent.  ARMCI matches
the GA model well and allows GA to support arbitrary dimensional
arrays efficiently.


\section{Implementing Adlib Communication in ter\-ms of
ARMCI\label{sec:AdlibOnARMCI}}

As explained in section \ref{sec:schedules}, the kernel Adlib library
is implemented directly on top of MPI.  We anticipated that
reimplementing parts of the collective communication library on top of
ARMCI would lead to improved performance on shared memory platforms (at
least).  Also, this task exercises the ARMCI API in a new context and
should help to prove it as a good interface.

There is a difference between the way the GA array descriptor is
implemented and the way distributed arrays are represented in Adlib,
arising from the different assumptions about communication: remote
memory access versus collective communication.  In Adlib the basic DAD
need only be accompanied by a single pointer to the base address
for the locally held elements.  In Global Arrays there is locally-held
{\em table of pointers} containing base addresses where the array
segments are stored in every peer process.

To exploit ARMCI, an Adlib program must at some point set up this table of
remote pointers.  In principle (at least in non-shared memory
implementations) this could be done either when the distributed array
was created, or at the time the remote pointers are needed (at schedule
construction or execution time).  The advantage of the latter approach
is that the old Adlib interface might be preserved---the user code
could maintain just the local base address, and remote address tables
could be regenerated on the fly.  But regenerating the table clearly
involves some rendezvous between the peers, which is likely to involve
a significant and repeated overhead.  Also it would imply some
complication of the ARMCI interface, because currently the {\tt
ARMCI\_Malloc} function is responsible for both allocating memory for
array elements and distributing handles to peers.  These two functions
would have to be separated.

These considerations indicate that the table of remote pointers should
be set up at array creation time (as in GA).  As explained in section
\ref{sec:DAD} the kernel DAD structure does not hold pointers to array
data, so this requirement does not imply any change to the
organization of that structure.  It {\em does} imply some changes at
the level of higher-level ``container classes'' (such as the ad++ {\tt
Array}$R${\tt<T>} classes) and---more importantly---in the interface to the
schedule classes.  The constructors (or {\tt execute} members) of the
schedules must be passed tables of pointers rather than individual
pointers.

\subsection{New classes vs reimplementation of existing
ones\label{sec:newClasses}}

As noted above, the ARMCI implementation implies changes to the actual
interfaces of kernel Adlib classes.  We might consider changing the
Adlib API in future releases, but this wasn't practical on the time
scale of the current project.  In any case it seems that there will
continue to be a place for arrays that {\em do not} carry full remote
table pointers.  Creating a standard Adlib distributed array is a
purely local operation and therefore very fast.  In the data parallel
style of programming it is quite common to create short-lived
distributed arrays as temporaries; to properly support this paradigm we
want creation of these arrays to be as fast as possible.  Creating a GA
array is itself a collective operation, and intrinsically slower.

In summary, we can distinguish two useful kinds of distributed
array: the ``local'' distributed array of Adlib and the ``global''
distributed array of GA.

For the initial experiments we added a new series of container classes to
ad++ that parallel the pre-existing ``local'' distributed
arrays, but implement GA style ``global'' distributed arrays:
\small
\begin{verbatim}
        template<class T>
          class ArrayRM1 ;

        template<class T>
          class ArrayRM2 ;

        template<class T>
          class ArrayRM3 ;

        ...
\end{verbatim}
\normalsize
(In fact these classes are derived from the corresponding classes in
standard ad++, so ordinary Adlib operations that don't need remote
pointer tables are still applicable to them.)

Similarly for communication schedules, we have implemented
a subset that parallels standard Adlib ones, with different
APIs.  For example
\small
\begin{verbatim}
  class RemapRM {
  public :
    RemapRM(const DAD* dst, void* dstDat,
            const DAD* src, const int srcOff, void** srcPtrs,
            const int len) ;

    void execute() ;
  } ;
\end{verbatim}
\normalsize
For further discussion see section \ref{sec:ARMCIBasedCollective}.

In the longer term we may want to consider a more uniform treatment of
the two kinds of distributed array.  One possibility is to elevate the
{\tt group} field of the DAD to include more general contextual
information, in a manner reminiscent of the way MPI decorates its
groups to form communicators.  This contextual information might
include information about whether the array supports remote memory
operations, the associated extra tables, and perhaps other flags
specifying information such as whether an array is held in-core or
out-of-core.

\subsection{ARMCI-compatible Adlib arrays}

One other minor problem was encountered in the construction
of Adlib-like arrays supporting ARMCI-based communication.

In the Adlib model it is supposed to be possible for {\em any} group of
processes to collectively create a distributed array.  The elements are
divided across processes in this group, which may or may not include all
processes executing the program as a whole.  Processes {\em outside}
the group do not have to engage in the array creation at all---they may
be executing some completely different code at the time.  This feature
of the library was originally introduced to support translation of a
kind of nested parallelism that can occur in HPF programs when (for
example) some intrinsically parallel Fortran 90 array expressions
appear nested inside a {\em forall} construct.  It has been used in
practise to support ``multiple data parallelism'' \cite{HPFKelp}, and
also in work on a parallel database classification algorithm using a
divide-and-conquer style of parallelism, combined with data parallelism
\cite{Classification}.

The ARMCI functions {\tt ARMCI\_Malloc}, {\tt ARMCI\_free} allocate raw
shared or remotely accessible memory.  In the original definition they
were collective across all processes on which the main program is
initiated.  To support the Adlib requirements we prototyped
versions with interface
\small
\begin{verbatim}
  ARMCI_MPI_Malloc(MPI_Comm comm, void *ptr_arr[], int bytes) ;
  ARMCI_MPI_Free(MPI_Comm comm, void *ptr) ;
\end{verbatim}
\normalsize
These versions are collective only within the process group associated
with the communicator argument.  Similar functions will now be
incorporated in the ``official'' ARMCI API.

To illustrate the how construction of an Adlib array is affected,
we sketch the ad++ constructors for the original {\tt Array2<T>}
and the new {\tt Array2RM<T>} class\footnote{We re-emphasize that
ad++ is just one possible, particularly simple, interface to kernel
Adlib.  It illustrates the issues.}.
The original constructor is
\small
\begin{verbatim}
  Array2(const Range& x_0, const Range& x_1, const Group& p) : 
                                   Section2<T>(p) {
    int size = 1 ;

    maps [0] = Map(x_0, size) ;
    size *= x_0.volume() ;
    maps [1] = Map(x_1, size) ;
    size *= x_1.volume() ;

    data = new T [size] ;
  }
\end{verbatim}
\normalsize
The total size of the local segment is accumulated in the {\tt size}
variable.  The {\tt volume} member of {\tt Range} returns the number
of locally held index values for each dimension.
Intermediate values of {\tt size} are memory strides associated 
with array dimensions and these are stored with the ranges in the
{\tt Map} objects of the DAD.  Finally the space for the local elements
is allocated.  The destructor just deletes {\tt data}.

The ARMCI-compatible class contains two
extra fields, {\tt ptrs}---the table of remote pointers---and
{\tt offset}---a constant offset to be added to these pointers (it may
become non-zero for certain array sections constructed from the
original array).  The constructor is
\small
\begin{verbatim}
  ArrayRM2(const Range& x_0, const Range& x_1, const Group& p) : 
                                     SectionRM2<T>(p) {
    int size = 1 ;

    maps [0] = Map(x_0, size) ;
    size *= x_0.volume() ;
    maps [1] = Map(x_1, size) ;
    size *= x_1.volume() ;

    ptrs = new void* [AdlibNProcs()] ;
    malloc_rm(ptrs, &comm, sizeof(T) * size) ;
    data = (T*) ptrs [AdlibLocId()] ;

    offset = 0 ;
  }
\end{verbatim}
\normalsize
The difference is that {\tt ptrs} must be allocated, and then there is
a call to {\tt malloc\_rm}.  This is a wrapper for {\tt ARMCI\_MPI\_Malloc}.
The {\tt malloc\_rm} call sets up an MPI communicator for the Adlib active
process group and passes it to {\tt ARMCI\_MPI\_Malloc}.  Note that the
{\tt comm} argument is actually an {\em out} argument; the communicator
is cached in the {\tt ArrayRM2} record so that it can eventually be
passed to a matching {\tt free\_rm} call in the \verb$~ArrayRM2$
destructor, and thus to {\tt ARMCI\_MPI\_Free}.  The destructor also
deletes {\tt ptrs}.

\subsection{ARMCI-based collective communication
schedules\label{sec:ARMCIBasedCollective}}

A representative selection of the Adlib communication
schedule classes were reimplemented in terms of ARMCI.  These included
the {\tt Remap} class, which implements copying between regular
sections of distributed arrays, the {\tt Gat\-her} class, which
implements copying of a whole array from a source array indirectly
subscripted by some other distributed arrays, and a few related
classes.

As explained in section \ref{sec:newClasses} the ARMCI-based schedules
needed a different API to the standard Adlib schedules.  Different
arguments are needed for constructor and {\tt execute} member.  In the
event we implemented a completely new set of low-level schedules based
on ARMCI.  They are listed in Table \ref{tab:ARMCIschedules}, which
should be compared with Table \ref{tab:schedules}.  Compared with the
MPI-based schedules these classes have a simple role (and
implementation).  They just store lists of individual data requests.
The execution members simply dispatch the ARMCI operations indicated
by these lists.  The only difference between the collective and
non-collective versions is that the execution members of the collective
versions incorporate barrier synchronizations and, where appropriate,
{\tt ARMCI\_AllFence} operations.

\begin{table}[tb]
\caption{Low-level ARMCI-based schedules\label{tab:ARMCIschedules}}
\vspace{0.1in}
\centerline{
\begin{tabular}{| c | c | c |}
\hline
                & operations on ``words'' & operations on ``blocks'' \\
\hline
non-collective  & {\tt RMSchedule }       & {\tt BlockRMSchedule } \\
\hline
collective      & {\tt CollRMSchedule }   & {\tt CollBlockRMSchedule } \\
\hline
\end{tabular}
}
\end{table}

The {\tt RMSchedule} and {\tt CollRMSchedule} schedules use the
ARMCI operations {\tt ARMCI\_GetV} or {\tt ARMCI\_PutV} in their
execution members.  The {\tt BlockRM\-Schedule} and {\tt
CollBlockRMSchedule} schedules generally use {\tt ARMCI\_GetS} or {\tt
ARMCI\_\-PutS} (recall that the ``blocks'' recognized by low-level Adlib
schedules are multidimensional local array sections, parametrized by
vectors of extents and memory strides).

New top-level schedules are implemented in terms of these
lower level ones.  The interface of {\tt RemapRM}, for example, was
already displayed in section \ref{sec:newClasses}:
\small
\begin{verbatim}
  class RemapRM {
  public :
    RemapRM(const DAD* dst, void* dstDat,
            const DAD* src, const int srcOff, void** srcPtrs,
            const int len) ;

    void execute() ;
  } ;
\end{verbatim}
\normalsize
Where the old remap operation was implemented on top of the
point to point schedule {\tt BlockMessSchedule}, the natural
basis of {\tt RemapRM} is {\tt CollBlockRM\-Schedule}.  The detailed
API and usage of {\tt CollBlockRMSchedule} is rather different
from {\tt BlockMessSchedule}, so the old implementation of {\tt Remap}
in terms of {\tt BlockMessSchedule} cannot easily be recycled.  Luckily
Adlib includes a generalization of {\tt remap} called {\tt vecGather}.
The schedule of the latter is implemented on top of the old ``remote
access'' schedule {\tt BlockDataSched\-ule}, whose usage 
is quite similar to {\tt CollBlockRMSchedule}.  So in practise {\tt
VecGather} was retargetted to {\tt CollBlockRMSchedule} to create an
ARMCI-based schedule called {\tt VecGatherRM}.  The resulting code was
simplified to create {\tt RemapRM}.

A minor difference between the the new {\tt RemapRM} schedule and
the old {\tt Remap} schedule (and conventional Adlib schedules
in general) is that source and destination data pointers are passed to 
the {\em constructor} rather than the execute member.  This is
because the lists stored in the ARMCI-based schedules are essentially
arguments for ARMCI functions, which include pointers.  Conventionally
Adlib stored only offsets (not absolute data addresses) in schedules.
This is less convenient when targetting ARMCI\footnote{A goal
in developing Adlib is always that it should be as ``lean'' as possible.
Of course we could have preserved more features of the old API.
They were changed to avoid superfluous array allocation and copying
inside the implementation.}.

As discussed earlier, the table of remote pointers for the source array
must be passed in to the schedule.  Only the local pointer for the
destination array is needed.  The argument {\tt srcOff} is the offset
mentioned in the previous section.  This may be non-zero if the source
is an array section.

\subsection{Benchmarking the ARMCI-based schedules}

The benchmark presented here is based on the {\em remap} operation.
The particular example chosen abstracts the communication in the
array assignment of the following HPF fragment
\small
\begin{verbatim}
          real a(n, n), b(n, n)
    !hpf$ distribute a(block, *) onto p
    !hpf$ distribute b(*, block) onto p
          a = b
\end{verbatim}
\normalsize
The destination array is distributed in its first dimension and the
source array is distributed in its second dimension, so the assignment
involves a data redistribution requiring an all-to-all communication.
In practice the benchmark was coded directly in ad++ (rather than
Fortran) and run on four processors of an SGI Power Challenge.  The
timings for old and new versions are given in Figure
\ref{fig:AdlibPerf} and table \ref{tab:AdlibPerf}.  Generally the
implementation of the operation on top of ARMCI is more than twice as
fast as the MPI implementation\footnote{For small arrays the MPI
version was actually slightly faster due to extra barrier
synchronizations used in the initial ARMCI implementation of {\tt
RemapRM}.  For now these synchronizations were implemented naively in
terms of message-passing, and contribute around 190$\mu$s to the total
times.}.  The MPI implementation is MPICH using the shared memory
device, so the underlying transport is the same in both cases.

\begin{figure}[tb]
\centerline{\psfig{figure=remap.eps,width=4in}}
\caption{Timings for original MPI vs new ARMCI implementation of {\em remap}.
The operation is a particular redistribution of an $N$ by $N$ array.
\label{fig:AdlibPerf}}
\end{figure}

\begin{table}[tb]
\caption{Times in milliseconds plotted in Figure \ref{fig:AdlibPerf}
\label{tab:AdlibPerf}}
\vspace{0.1in}
\centerline{
\small
\begin{tabular}{lll}
Array extent &  MPI-based  & ARMCI-based \\
\\
32           &  0.20       &  0.21       \\
64           &  0.43       &  0.25       \\
128          &  1.11       &  0.52       \\
256          &  2.07       &  1.08       \\
512          &  11.7       &  5.35       \\
768          &  29.4       &  13.8       \\
1024         &  62.4       &  27.2       \\
1280         &  86.8       &  39.0       \\
1536         &  127.8      &  56.1       \\
\end{tabular}
}
\end{table}



\subsection{One-sided communication schedules for Adlib\label{oneSidedAdlib}}

By replacing the {\tt CollBlockRMSchedule} base class of
{\tt RemapRM} with the non-collective version {\tt BlockRMSchedule}
(which doesn't include barrier synchronizations), and
stripping down the code so that the destination array is
just a local sequential array rather than a distributed
array, we produced a schedule called {\tt Get}.  This
schedule supports one-sided ``get'' (also ``put'') operations
functionally equivalent to the GA operations {\tt ga\_get} and
{\tt ga\_put}.

The ad++ {\tt get} functions based on the {\tt Get} schedule have
interfaces like:
\small
\begin{verbatim}
  template<class T>
  void get(const int str_0, const int str_1, T* base,
           const SectionRM2<T>& src) ;
\end{verbatim}
\normalsize
The destination local array is parametrized by a base address and a
set of local memory strides.  The last argument is an arbitrary regular
section of a global distributed array.  In particular this can be an
arbitrary rectangular patch, as for {\tt ga\_get}.

As a test case we transcribed the {\tt jacobi.F} example in the
GA release to C++ using ad++ {\tt ArrayRM2} arrays and the {\tt get}
operation, and verified that it reproduced the same results as the
original Fortran.

In the past the lack of this kind of one-sided communication has been a
notable gap in the functionality of Adlib.  ARMCI will allow us to
close this gap in future releases.


\section{Adapting Global Arrays to use the Adlib array
descriptor\label{sec:GAonDAD}}

In an orthogonal but complementary development, we have produced a PCRC
version of GA that internally uses a PCRC/Adlib array descriptor to
maintain the global array distribution parameters.  The modularity of
the GA implementation, enhanced by the introduction of the ARMCI layer,
made changing the internal array descriptor a relatively
straightforward task.  An immediate benefit is that direct calls to the
optimized Adlib collective communication library will be possible from
the modified GA.

\subsection{Extensions to the DAD}

The original Adlib array descriptor supported all distribution formats
specified in HPF 1.0 \cite{HPFStandard}.  These formats include simple
block, simple cyclic, and block-cyclic distributions.  It did not
support the kind of irregular block distributions allowed in Global
Arrays.  Version 2.0 of HPF \cite{HPF96} specified equivalent block
irregular distributions as an ``approved extension'', so provisions had
been made in Adlib for this extension.  But the extension had not been
implemented at the time the current project started.

A new subclass of {\tt Range} called {\tt IrregRange} was added
to the hierarchy of figure \ref{fig:ranges}.  The documented interface
is
\small
\begin{verbatim}
  class IrregRange : public Range {
  public :
    IrregRange(Dimension dim, int blocks []) ;
  } ;
\end{verbatim}
\normalsize
The constructor creates an irregular block-distributed range
distributed over the process grid dimension {\tt dim}.  The array {\tt
blocks} has $P$ elements where $P$ is the extent of the grid
dimension.  These elements define the block size associated with each
process.  The extent of the range is
\small
\begin{displaymath}
\sum_{i = 0}^P \mbox{\tt blocks [} i \mbox{\tt ]}
\end{displaymath}
\normalsize
This interface for the constructor is a natural one for directly
supporting the new {\tt GEN\_BLOCK} distribution format in HPF 2.0.
We also included a second (undocumented) constructor that is slightly
more convenient for implementing Global Arrays API:
\small
\begin{verbatim}
    IrregRange(const int extent, Dimension dim, Integer map []) ;
\end{verbatim}
\normalsize
The definition of {\tt map} agrees with the definition of the
corresponding arguments of {\tt GA\_CREATE\_IRREG}: the elements are
the smallest global index held in each process.

Internally {\tt IrregRange} uses an algorithm borrowed from the
GA implementation for address translation.  In particular the
inline function {\tt find\-Block}, used in the global-to-local
subscript conversion method called {\tt locat\-ion}, is an adapted
version of the corresponding GA macro.

There remains one significant difference between the layout of elements
in a GA array and the corresponding Adlib array, even using {\tt
IrregRange}.  In Adlib the same volume of memory is allocated in every
member of the process group specified for the array:  this volume is
the product of the {\em maximum} local block sizes from each range.  GA
allocates a volume equal to the product of the local block sizes ({\em not}
the maximum), and this volume is usually different in each
process.  Adlib exploits the constancy of the volume, assuming that
memory strides for remote array segments are the same as for the local
segment.  These strides are stored in the DAD.  Rather than try to modify
the DAD it was decided to live with this difference.  Although it is
quite likely that some Global Arrays programs rely on the assumption
that the locally allocated array volume is identical to the number of
locally-held global array elements, we will see later that we did not
encounter problems with the example programs in the GA release
package.


\subsection{Modifications to Global Arrays}

The original plan called for implementation of a small subset of
Global Arrays, including {\tt ga\_create}, {\tt ga\_destroy},
{\tt ga\_get}, {\tt ga\_put} and {\tt ga\_acumulate}.
Originally it was intended that these would be Fortran wrappers around
calls to the {\tt Get} schedule described in section 
\ref{oneSidedAdlib}.

To validate the subset, one or more test programs would be needed.
Ideally these would be semi-realistic Global Arrays programs.
All the examples in the GA release use around a dozen different GA
operations.  Unless some very minimal test code was invented,
it seemed that the proposed ``small subset'' would have to be 
bigger than the original plan assumed.  It became clear
that it would be very useful to recycle code from the existing Global
Arrays implementation to create the new interface.  Most of the
relevant routines are implemented in the original GA source file {\tt
global.core.c}.
That file was about 3000 lines long.  It was hoped that by
cutting down the number of entry points to just those needed for the
current study, and retargetting to Adlib with ARMCI/MPI, a much smaller
file would be produced.  Before that attempt got very far, PNNL
provided a new implementation of core global arrays based
on ARMCI.  The corresponding file is called {\tt global.armci.c}.

The plan now altered slightly.  Rather than providing wrappers around
the new Adlib {\tt Get} schedule, we would modify
the file {\tt global.armci.c} to use the Adlib DAD, but leave direct
calls to ARMCI in the file.
The principle functions that would
need modification were the ones listed in table \ref{tab:modified}.
The structure {\tt global\_array\_t} from GA would be retained, because
it incorporates various fields used for management of arrays that
have no counterpart in the Adlib DAD.  But several fields---{\tt dims},
{\tt chunk}, {\tt nblock}, {\tt scale}, {\tt mapc}, {\tt lo}---would
be deleted and replaced by a DAD.

\begin{table}[tb]
\caption{Functions from {\tt global.armci.c} requiring
modification\label{tab:modified}}
\vspace{0.1in}
\centerline{
\small
\tt
\begin{tabular}{|l|l|l|}
\hline
nga\_create         & nga\_locate\_          & ga\_scatter\_acc\_local \\
nga\_create\_irreg  & nga\_locate\_region\_  & ga\_scatter\_           \\
ga\_create\_irreg   &                        & ga\_scatter\_acc\_      \\
ga\_duplicate       & ga\_fill\_             & ga\_gather\_local       \\
                    &                        & ga\_gather\_            \\
nga\_distribution\_ & nga\_proc\_topology\_  &                         \\
ga\_access\_        &                        & nga\_read\_inc\_        \\
                    & nga\_put\_             &                         \\
                    & nga\_get\_             &                         \\
                    & nga\_acc\_             &                         \\
\hline
\end{tabular}
}
\end{table}

Various kinds of changes were needed to {\tt global.armci.c}
and the associated header {\tt global.armci.h}.  A number
of global changes relate to need to compile these files as C++
rather than C.  Many of the functions had old-style C argument
lists, and these had to be changed to ANSI-style.  Function
prototypes had to be added in various place, or the function
definitions had to be reordered to avoid forward references
without prototypes.  In some places {\tt extern "C"} clauses
were needed.  In some places explicit casts were needed to
pass the more stringent type-checking of C++.

Various old functions and macros became irrelevant and were
deleted, including {\tt findblock}, {\tt gam\_GetRangeFromMap}, {\tt
gam\_CountElems}, {\tt gam\_Compute\-PatchIndex}, {\tt gam\_ComputeCount}
and {\tt gam\_setstride}.

A few new classes and functions were added.  A class similar to
the {\tt Get} schedule, but specialized to reproduce exactly the
functionality of {\tt nga\_get}, etc, was defined (lists of block
moves had to be constructed anyway, because of the way GA randomizes
the order of of individual remote accesses).  The class is called {\tt
PatchTransfers}.

Finally, functions in table \ref{tab:modified} had to
be modified to create and use the DAD fields instead of the
deleted {\tt global\_array\_t} fields.  In {\tt ga\_create\_irreg} 
a suitable Adlib process grid has to be constructed or found before
setting up the DAD.  The initialization and finalization functions had
to initialize and finalize Adlib.

Various other header files included by {\tt global.armci.c}
also had to be modified to meet the strict prototyping
requirements of C++.  These included 
{\tt globalp.h}, {\tt global.h} and {\tt ma/macdecls.h}.

While this work was in progress PNNL provided a beta version of the new
n-dimensional GA.  It was updated with the NPAC modifications, and our work
continued from that basis.

Once the new version of the core GA file could be compiled, we had
to address the problem of linking it.  One option was to try to extract
a suitable subset of necessary modules in the GA release (including,
for example, the MA memory allocator).  Selecting a suitable subset
seemed more difficult than recompiling the whole GA release to use the
new core module.  So the make files were modified to compile
just the core module with C++ and link it with the rest of the
library.

After some debugging, all relevant examples in the {\tt global/testing}
dir\-ectory---
{\tt test.x},
{\tt jacobi.x},
{\tt patch.x},
{\tt perf.x},
{\tt perfmod.x},
{\tt perform.x},
{\tt p.x},
{\tt t.x},
{\tt testeig.x},
{\tt testsolve.x} and
{\tt ndim.x}---could be run successfully with the modified GA.
A more detailed list of changes to the standard version of GA is available
from NPAC.

\subsection{Benchmarking the DAD-based GA primitives}

Figures \ref{fig:locGet2d} and \ref{fig:remGet2d} and tables
\ref{tab:locGet2d} and \ref{tab:remGet2d} give timing results
for the modified version of {\tt ga\_get}, compared with
the (new ARMCI-based) version of GA provided by PNNL.
Timings were obtained using the program {\tt global/testing/perf.F}
from the GA release (once again on an SGI Power Challenge).

We conclude that using the object-oriented DAD of Adlib introduces a
small overhead in the {\tt ga\_get} operation.  This overhead appears to be
generally less than a microsecond%
\footnote{The first naive implementation of GA using the DAD introduced
a larger overhead (a few microseconds).  Some optimizations were applied
to reach the performance quoted here.  However the current version is
still not really aggressively optimized---the work stopped when
overall performance was comparable with the original GA.}%
.

\begin{figure}[tb]
\centerline{\psfig{figure=locGet2d.eps,width=4in}}
\caption{Timings for modified Global Arrays vs original.
The operation is a {\em get} operation on a $dim$ by $dim$ patch of an
array.  Timings when patch is locally resident.
\label{fig:locGet2d}}
\end{figure}

\begin{table}[tb]
\caption{Times in microseconds plotted in Figure \ref{fig:locGet2d}
\label{tab:locGet2d}}
\vspace{0.1in}
\centerline{
\small
\begin{tabular}{lll}
Patch width  &  Using DAD  & Original GA \\
\\
1            & 5.69        & 5.23        \\
3            & 7.20        & 6.22        \\
4            & 6.45        & 6.03        \\
9            & 9.66        & 8.33        \\
16           & 13.9        & 13.3        \\
24           & 22.8        & 21.5        \\
30           & 33.0        & 29.1        \\
48           & 70.7        & 69.1        \\
64           & 138         & 126         \\
91           & 351         & 337         \\
128          & 835         & 825         \\
171          & 2021        & 1967        \\
256          & 7513        & 7433        \\
353          & 15930       & 15766       \\
440          & 25780       & 25638       \\
512          & 26623       & 34371       \\
\end{tabular}
}
\end{table}

\begin{figure}[tb]
\centerline{\psfig{figure=remGet2d.eps,width=4in}}
\caption{Timings as for Figure \ref{fig:locGet2d},
but in case where patch is held remotely.
\label{fig:remGet2d}}
\end{figure}

\begin{table}[tb]
\caption{Times in microseconds plotted in Figure \ref{fig:remGet2d}
\label{tab:remGet2d}}
\vspace{0.1in}
\centerline{
\small
\begin{tabular}{lll}
Patch width  &  Using DAD  & Original GA \\
\\
1            & 11.6        & 10.5        \\
3            & 14.2        & 13.7        \\
4            & 15.3        & 15.0        \\
9            & 21.8        & 20.7        \\
16           & 63.8        & 63.0        \\
24           & 111         & 114         \\
30           & 144         & 143         \\
48           & 373         & 368         \\
64           & 636         & 640         \\
91           & 1276        & 1247        \\
128          & 2370        & 2075        \\
171          & 4215        & 4449        \\
256          & 12936       & 11970       \\
353          & 23314       & 23443       \\
440          & 34540       & 33341       \\
512          & 39244       & 40716       \\
\end{tabular}
}
\end{table}

\clearpage

\bibliographystyle{plain}
\bibliography{report}

\end{document}


