EXPLORING  THE  POTENTIAL  OF  THE  SMP  CLUSTER  ARCHITECTURE
ON  THE  SILICON  GRAPHICS  POWER  CHALLENGE  ARRAY
Paul R. Woodward
University of Minnesota
Director, Laboratory for Computational Science & Engineering
August 9, 1996
I.   Introduction
This is a report of activities in the Laboratory for Computational Science & Engineering (LCSE) in collaboration with industrial sponsors and government sponsored laboratories.  The industrial collaborators are Silicon Graphics (SGI), Seagate Technology, Ciprico, Micro Technology Inc. (MTI), Ampex Data Systems, DEC Networking Division, Prisa, and GeneSys.  The government sponsored laboratories are the National Center for Supercomputing Applications (NCSA), Los Alamos National Laboratory (LANL), the Army Research Laboratory (ARL), Lawrence Livermore National Laboratory (LLNL), the Electronic Visualization Laboratory (EVL), the Minnesota Supercomputer Institute (MSI), the Pittsburgh Supercomputer Center (PSC), NASA Goddard Space Flight Center, the Jet Propulsion Laboratory (JPL), the Cornell Theory Center, Argonne National Laboratory (ANL), and the White Sands Missile Range.  The industrial participation has generally involved loans or donations of equipment, often preproduction equipment, provision of prerelease software, and direct involvement of personnel.  Results of these collaborations are often demonstrated at exhibits around the country as well as at the LCSE.  Participation from government sponsored laboratories most commonly involves direct involvement of personnel, although grants and contracts, whether from a parent federal agency or directly from the laboratory involved, are common additional mechanisms.  Occasionally, as with JPL and White Sands, travel and other incidental project expenses are provided.
These collaborative projects are described in this single document because they share a common theme, scientific computation and visualization using the SMP cluster architecture, the Silicon Graphics Reality Engines, the very high resolution PowerWall display, and support of these systems with high speed RAID disk and tape systems.  The LCSE is focussed upon developing computational techniques on SMP clusters (and distributed shared memory, or DSM, machines, as they emerge) and on developing the infrastructure to explore very large data sets interactively which are produced by supercomputer computations on this architecture.  The LCSE is a development lab;  it is not a production facility.  The primary output of the lab consists of new computational techniques, new computational, visualization, and data storage systems and system software, and demonstrations of the effectiveness of these developments through their creative use.  Formal projects, grants, and agreements which make this work possible are listed at the end of this document.  The account which follows attempts to present these activities in a more unified format.
II.   Equipment Configuration
During the last year the equipment configuration in the LCSE has changed considerably as a result of equipment loans from industrial partners.  For about 2 months in advance of the Supercomputing ’95 conference and exhibit in San Diego last December, a special equipment configuration allowed us to develop demonstrations for both the conference-organized events and for the Silicon Graphics booth.  This special equipment configuration is shown in Figure 1.  Its central feature was the cluster of 4 powerful SGI machines, an 8-processor Power Onyx on 9-month loan to the LCSE since its founding, a 12-processor Power Challenge machine provided by the Army Research Laboratory to the LCSE for a two-year project in distributed virtual memory research, and two 8-processor Power Onyxes loaned to the LCSE by SGI to drive our PowerWall display in preparation for its role in the SGI booth at Supercomputing ’95.  For the Supercomputing ’95 event we worked with LCSE industrial partners Seagate Technology and MTI to construct a terabyte disk subsystem to drive the PowerWalls, CAVE, and Immersadesk systems at the show.  This terabyte disk system, the first to our knowledge to be set up as a single file system on a single machine, our Power Challenge machine in the LCSE, was the centerpiece of a technology open house which we held in the LCSE the second week after the Supercomputing ’95 show.
Figure 2 shows the present equipment configuration in the LCSE, including some items which have either been ordered or which are likely to be included in our SGI equipment loan for this year.  Through an SGI equipment loan, our 12-processor Power Challenge machine has been augmented to include the maximum 18 processors.  We hope to add 16 R-10000 processors to the Onyx machine in this diagram during the next year, so that we will once again have a “cluster” of at least two machines in the lab.  The Onyx has two Infinite Reality graphics engines, each with two raster managers.  This single chassis drives the four panels of the 3200x2400 pixel PowerWall display.  This figure indicates new fibre channel disks and disk arrays which we have either ordered or budgeted for the coming year.  The Ciprico fibre channel RAID systems will provide the principal high speed storage for the PowerWall display, and we intend to use them also for experimentation with out-of-core computation strategies.
III.   The Fortran-P Project:  Template Codes and the PPMLIB Library
This project is funded by the DoE, the ONR, and to a limited extent through the LCSE MetaCenter Regional Alliance grant from the NSF.  This project involves the development of tools to aid in the writing of application codes for massively parallel computers (MPPs) and for clusters of shared memory multiprocessors (SMPs), like the Power Challenge Array.  This project grew out of the very successful code parallelization experiments carried out by my group and that of Matt O’Keefe at the AHPCRC.  The project now finds its home in the LCSE, as is appropriate given our lab’s focus on the new computer architecture represented by SMP clusters like the Power Challenge Array.
The Fortran-P project originally focused upon the generation of automatic code translation tools, the Fortran-P precompiler, to parallelize a stylized Fortran-77 source code for execution on MPPs like the SIMD Connection Machine CM-5 of the AHPCRC.  Due to the strong limitations of the SIMD programming model (the data parallel model) which we were unable to overcome on the CM-5, the original Fortran-P precompiler for that machine is an imperfect starting point for a more general and far more powerful tool targeted at SMP cluster programming.  Although at the outset of this project we had many ideas for the extension of the data parallel Fortran-P precompiler to address the much broader class of applications which we believe can be efficiently implemented on the SMP cluster architecture, we decided that it was best to work through a variety of specific examples before fixing this expanded Fortran-P syntax and semantics.  We therefore devoted considerable energy during the last year to the development of template application codes.  These can serve as well documented and cleanly implemented examples to others wishing to port applications to systems like the SGI Power Challenge Array.  These template codes have now been fairly well documented in the articles appended to this report, and we hope that they will now prove more useful as a result.  We have also used these templates and the experience of Matt O’Keefe’s group working with the more complex MICOM ocean circulation code to determine a draft specification of Fortran‑P.
We have developed 4 template implementations of our PPM gas dynamics code.  Two of these four template codes are 2-D and two are 3-D codes.  We have tried to design them in such a way that their structure and the multitasking and network parallelization techniques which they illustrate are much more general than the specific gas dynamical applications which they embody.  We have shared these template codes with collaborators in areas as diverse as aeroelasticity, ocean modeling, mesoscale weather prediction, homogeneous turbulence, laser fusion, cosmological gas dynamics, solar convection modeling, and general relativity.  We are also beginning a new collaboration with Los Alamos to investigate the use of the techniques in this code template for radiation transport calculations.  Many of these collaborations fall under funded projects in our research group or loose associations with scientists at DoE and NASA laboratories.  We have also offered this template technology to our colleagues at ARL.  The template codes which we have developed to date are described briefly below.  Some of them can be found on the LCSE Web site at www.lcse.umn.edu.  Detailed documentation associated with these templates is appended to this report.
IIIa.   The Tiled 2-D sPPM Template Code
This is our most popular and perhaps most general template code.  It represents an implementation of the simplified PPM gas dynamics algorithm, sPPM, which attempts to address the issues involved with porting data parallel programs developed for machines like the CM-5 Connection Machine, the Cray T3D, or the IBM SP-2 to the SMP cluster architecture as represented by the SGI Power Challenge Array.  This template was developed expressly for demonstration at the October,1995, LCSE workshop on SMP cluster programming.  In the process of inviting several of the participants to that workshop, it came to my attention that this particular porting problem was shared by many investigators, and few had either the experience with the new architecture or the time to develop effective porting strategies.  This template code illustrates several essential points in dealing with the SMP cluster architecture, and feedback which we have received from our collaborators who have been working with this template code has assured us that this template serves indeed as a useful guide to the code porting process.
For this tiled 2-D template code, we built a special version of our sPPM gas dynamics module in order to simulate the structure of many other codes in use in the HPC community.  We built a version of sPPM which works directly on an Eulerian grid, and we packaged it in such a way that a fundamental unit of computation, invoked by a single subroutine call, consists of updating a rectangular tile of a larger 2-D grid.  Although this code is 2-D for simplicity, its structure is similar to many 3-D codes in the climate and meteorological community, where problem domains are typically very thin in one of the 3 dimensions.  This template is therefore guiding the efforts of collaborators at the University of Colorado, NCSA, the University of Oklahoma, and the University of Miami in porting ocean and atmosphere simulation codes to SMP clusters.  This last effort is in collaboration with Matt O’Keefe and Aaron Sawdey at the LCSE.  Their experience in applying this template to the MICOM ocean code has led to Aaron Sawdey’s Ph.D. thesis topic proposal, the basis of his successful preliminary oral exam this spring.  Aaron will be working on the Fortran-P translator project.
This tiled 2-D sPPM template code incorporates many code transformations within the sPPM module which we have found over the years to yield outstanding floating point performance on RISC microprocessors.  In fact, we have clocked the performance of this sPPM module at 131 Mflop/s (single precision) and at 111 Mflop/s (double precision) on a single R-8000 CPU running at 90 MHz.  We used a run on the Cray C-90 of the Minnesota Supercomputer Institute in order to count the flops in this module.  The C-90 performance was 416 Mflop/s on this module.  The R-8000 performance was therefore roughly a quarter (double precision) to a third (single precision) of the C-90 performance.  This version of the sPPM algorithm involves 849 flops per cell update per 1-D pass, or 1698 flops per cell update in 2-D.  
This 2-D template code uses sophisticated multitasking techniques to remove unnecessary barriers and to make large sections of the computation entirely self-scheduled.  Communication between tasks within a single SMP is via shared memory, with no explicit message passing of any type.  Communication between SMPs over the HiPPI network is via the MPI message passing library.  This message passing is entirely asynchronous and comes at extremely low cost.  The overlapping of communication with computation in this template code is so aggressive that we used it at Supercomputing ‘95 to test our concept for transcontinental distributed computing (see below).  Performance figures for this tiled 2-D template code are quite good, and in fact represented both the best speed-ups and the best speeds reported for SMP cluster codes (on any manufacturer’s equipment) at our LCSE workshop, where competing work from NCSA and NASA Goddard, in particular, was reported.
The following table reports tests with the tiled 2-D sPPM template code carried out in the LCSE using 3 Power Onyx machines, each with 8 processors running at 90 MHz, with 2 GB memory, and using only a single HiPPI interface per machine for message passing.  The communication was carried out using our own, high performance version of a restricted subset of MPI.  The run times reported in the table below are for 100 time steps.  Speed-ups are relative to performing the run on a single processor, using the 768x256 grid and only a single tile.  Single processor times for the larger grids were obtained by scaling the results of this run, assuming the same efficiency in updating each grid cell.  The two runs on the smaller grids divided each SMP’s subdomain into a 4x4 set of grid tiles, while the two runs on the larger grids used 8x8 grid tiles for each SMP.  The unusually small performance increase in going from the 1536x512 to the 3072x1024 grid comes from changing the number of subtiles used on each SMP.  In both these cases, the efficiencies are limited principally by the amount of redundant computation which must be performed along the overlapping boundaries of the grid tiles.  The performance in Gflop/s reported in the table is calculated from the knowledge that 1698 flops are required to update each cell of the grid in 2-D.  No redundant or otherwise unnecessary flops are counted.  The Gflop/s figures in the table reflect the fact that performance of the sPPM gas dynamics module is degraded from its 1-D level when it is embedded in a 2-D application, where, for example, extra copying of the data and, worse, transposition of the data is required along with the fundamental gas dynamics computation.
Performance of the Tiled 2-D sPPM Template Code
on the SGI Power Challenge Array with 3x8 MIPS R-8000 CPU’s @ 90 MHz
GRID	time  (sec)	Speed-Up	Efficiency	Gflop/s	
768x256	21.04	17.6	73.3%	1.59	
1536x512	67.8	20.9	87.1%	1.97	
3072x1024	247.8	21.4	89.0%	2.16	
4608x1536	543.0	23.4	97.4%	2.21	

The following table reports results of tests using the tiled 2-D sPPM template code running in single-precision mode on the 6-machine array in the SGI booth at Supercomputing ’95.  These results were obtained under non-ideal conditions, and therefore do not represent a best effort.  Clearly, the grids chosen are far too small to expect good efficiency on this very powerful, 96-processor system.  For these tests, the MPI implementation MPICH, from Argonne National Laboratory, was used for the communications.

Performance of the Tiled 2-D sPPM Template Code
on the SGI Power Challenge Array with 6x16 MIPS R-8000 CPU’s @ 90 MHz
GRID	time steps	subtiling	NSMPx	NSMPy	Nthreads	time (sec)	
2304x768	250	8x8	3	2	16	96.6	
2304x768	250	8x8	3	2	8	167.9	
2304x768	250	12x12	3	2	16	101.6	
2280x760	250	10x10	3	2	16	91.6	
2304x768	250	8x8	1	1	8	847.4	

The speed-ups, efficiencies, and net performance figures corresponding to these test runs are given in the table below.  Speed-ups are once again relative to single CPU performance on a single tile of 768x256 cells, scaled to the relevant grid size.  For the 2304x768 grid, this scaled run time is, for 250 time steps, 8332 sec.  This indicates an efficiency of 123% for 8 CPU’s running this problem on a single SMP using 64 grid tiles.  This superlinear speed-up is presumably the result of increased cache efficiency when the grid is tiled.  It would be better to compare each multiprocessed run to the time required for a single CPU to perform the run in the most efficient manner, which could involve updating the grid tile by tile rather than all at once.  We have not had sufficient time to sort out all these subtleties, but we feel that our results are useful nonetheless.
Performance of the Tiled 2-D sPPM Template Code
on the SGI Power Challenge Array with 6x16 MIPS R-8000 CPU’s @ 90 MHz
GRID	subtiling	NSMPx X NSMPy	Nthreads	Speed-Up	Efficiency	Gflop/s	
2304x768	8x8	3x2	16	86.3	89.8%	7.78	
2304x768	8x8	3x2	8	49.6	103.4%	4.47	
2304x768	12x12	3x2	16	82.0	85.4%	7.39	
2280x760	10x10	3x2	16	89.1	92.8%	8.03	
2304x768	8x8	1x1	8	9.83	123%	0.886	

These results clearly indicate that it is indeed possible to extract 8 Gflop/s out of 96 R-8000 processors on a practical problem coded in a manner which can be emulated by a rather broad variety of other application codes.  We hope that this and our other template codes can help users of Power Challenge Array systems to realize the potential performance that they offer without having to go through the long development process which our team at the LCSE at the University of Minnesota is engaged in.  During the last 6 months, we used this template code for a test of distributed computing over the Internet which is described in a later section.  We also used it to obtain parallel speed-ups for the SGI machine from ARL augmented by 90 MHz processors loaned from SGI.  These speed-ups are shown in the figure below:

Speed-ups obtained with the 2-D sPPM Template Code on the 18-processor Power Challenge machine using a varying number of processors and size of grid.
8 of the CPUs run at 90 MHz, while the remaining 10 run at 75 Mhz
Speed-ups are reported relative to a 90 MHz processor.
The maximum possible speed-up is therefore 16.33

In this case, the speedups do not tell the whole story since they are relative to running the code on a single 90 MHz processor.  Since the system was equipped with a mix of 10 processors running at 75 MHz and 8 processors running at 90 MHz, the “maximum obtainable speedup” when running with N processors is:

	N						when 1 < N  8
	((8  90) + (N - 8)  75) / (N  90)  N	when 9  N  18

Therefore, a more telling statistic is the efficiency defined in this case to be the measured speedup divided by the maximum obtainable speedup.  These efficiencies are shown below:

2-D sPPM template code efficiencies measured when varying the number of processors used and the grid size.

These efficiencies, which are all above 95%, clearly show that the memory bus of the Power Challenge is by no means saturated by 18 R-8000 processors on a real application code if that code is properly structured.  It should also be noted that the self-scheduled construction of the calculation here has allowed a heterogeneous collection of processors to work together at high efficiency on this problem without special load balancing modifications to the code.
IIIb.   The Directionally Split 2-D sPPM Template Code
This template code is intended to represent a nearly optimally efficient implementation of a directionally split algorithm like PPM on the SMP cluster architecture.  As with the tiled 2-D template discussed above, this template incorporates a gas dynamical sPPM kernel which is highly optimized for the MIPS R‑8000 CPU, its large cache memory, and its relatively low bandwidth (160 MB/s per CPU) main memory.  In contrast to the tiled code just discussed, this code involves a minimum of redundant computation along subdomain boundaries.  This computational expense is incurred only at the boundaries along which the grid is subdivided among the SMPs in the cluster.  There is no more than a single internal boundary between grid tiles within the subdomain of each individual machine.  This feature makes the computation very efficient, even on rather small grids.
In this template code, as in the tiled code discussed above, communication is overlapped with computation in a very aggressive fashion.  In this template code, the overlapping is achieved by allowing a y-sweep of the computation to begin before the x-sweep has been completed.  The multitasking on each SMP is completely self-scheduled, with barrier synchronizations only after every two full time steps.  Messages are collected for dispatch in large blocks to neighboring SMPs on the network, and they are dispatched well before their receipt is required to continue the computation.  This technique allows enormous message passing latencies and very low message passing bandwidths to be tolerated without degrading overall performance.  This code has been thoroughly tested on a single cluster of SMPs, and it has also been used in tests of transcontinental distributed computing which are described below.  It should be noted that the sPPM kernel contained in this 2-D template breaks each 1-D sweep into separate Lagrangian hydrodynamics and remap computations.  This kernel runs at somewhat lower Mflop/s performance on the MIPS R-8000 than the sPPM kernel of the tiled 2-D template code.  Nevertheless, time to solution is not degraded, because the separation of the computation into Lagrangian and remap steps allows a larger time step to be used.
We have recently tested the performance of this template code on the ARL Power Challenge machine at the LCSE which has been augmented with additional R-8000 processors loaned to us by SGI.  The speed-ups thus obtained are given in the table below.  The times given in the table correspond to a single pair of time steps.  These results show a 5% initial drop in per-processor performance when the second processor is added to the job, but much smaller performance decreases as still more processors are added.  Considering the aggregate computing power of the 18 CPUs, these 2-D problems are quite small, and therefore the speed-ups reported in the table are quite acceptable.  As we shall see, the 3-D version of this template code achieves 97% efficiency with these 18 processors on a problem of reasonable size for this computing engine.
Performance of the Directionally Split 2-D sPPM Template Code
on the SGI Power Challenge Machine with 18 MIPS R-8000 CPUs
( 8 CPUs @ 90 MHz and 10 CPUs @ 75 MHz )
GRID	time  (sec)	CPUs	Speedup	Efficiency	
1024x1024	60.8	1	1	100%	
	32.2	2	1.88	94%	
	16.5	4	3.68	92%	
	8.50	8	7.15	89%	
	4.86	16	12.5	85%	
	4.47	18	13.6	83%	

Performance of the Directionally Split 2-D sPPM Template Code
on the SGI Power Challenge Machine with 18 MIPS R-8000 CPUs
( 8 CPUs @ 90 MHz and 10 CPUs @ 75 MHz )
GRID	time  (sec)	CPUs	Speedup	Efficiency	
2048x2048	249	1	1	100%	
	131	2	1.90	95%	
	67.0	4	3.74	93%	
	34.2	8	7.28	91%	
	19.5	16	12.8	87%	
	17.8	18	14.0	86%	

IIIc.   The Synchronous Communication 3-D sPPM Template Code
This 3-D sPPM template code was constructed for use by the DoE’s ASCI program (Accelerated Strategic Computing Initiative).  The ASCI teams at Los Alamos and Livermore required an illustrative example of an efficient SMP cluster implementation of a “real” code for use in their “ASCI Blue” project.  The “blue” component of the ASCI platform project focuses upon the SMP cluster architecture.  Since the SGI Power Challenge Array is the first representative of this new architecture in a commercial product, and because my LCSE team is well known for its pioneering work with this architecture, we were asked to provide the DoE with a template code based upon our PPM gas dynamics code.
We decided to keep the PPM gas dynamics module embedded in this 3-D template code as simple as possible.  However, we wanted the computation to memory access ratio for the simplified algorithm to reflect that of the full PPM code.  The simplified sPPM module in this benchmark is a compromise.  It omits several very complex computations contained in the PPM code while making up for their loss in part by more computational labor in its nonlinear Riemann solver.  The final result is a scheme which requires 817 flops (as measured by the Cray C-90’s flop counter) to update each cell of the grid for each 1-D sweep of the calculation.  This code is extensively commented, and it is based upon part of a course that I teach at the University, which is supported by rather detailed discussions in the course notes motivating the algorithm.  The implementation of this sPPM algorithm, which involves separate Lagrangian and remap steps, is optimized for cache-based microprocessors with software pipelining capabilities, like the MIPS R-8000 and several other recent chips or chip sets.  The sPPM module code executes at 442 Mflop/s on a Cray C‑90 processor.  (The identical code given to the R-8000 microprocessor executes more slowly on the Cray, for a number of reasons.  This less efficient Cray code executes at only 344 Mflop/s.)  In 32-bit mode, our machine from ARL, with 75 MHz R-8000 processors, produces 93.0 Mflop/s performance on this sPPM module code (the 90 MHz performance is 112 Mflop/s).  With 64-bit arithmetic, the R-8000 performance at 75 MHz is 82.9 Mflop/s.

Performance of the 3-D sPPM  ASCI  Template Code
on a 256x256x256 grid
on the SGI Power Challenge Machine with 12 MIPS R-8000 CPU’s @ 75 MHz
#  CPU’s	Speed-Up	Efficiency	Gflop/s	
1	1	100%	0.0921	
2	1.99	99.6%	0.183	
4	3.95	98.7%	0.364	
8	7.60	95.0%	0.700	
12	11.1	92.4%	1.021	

The implementation of this 3-D template code on a single SMP is aggressive.  The fundamental unit of work for a single CPU is the updating of a pencil of 8x8xN grid cells, where N is the length of the grid for this SMP’s subdomain in either the x, the y, or the z direction.  Special data layout and handling is incorporated in this template code to avoid any unnecessary data transpositions.  The tasks are completely self-scheduled, with barrier synchronizations only at the end of each 1-D sweep over the entire grid.  The multitasking efficiency can be seen from the above table of performance data generated in the summer of 1995 and reported at several conferences.
At the specific request of the DoE., we constructed this template code for synchronous message passing.  That is, at the end of each 1-D sweep over the grid, a barrier synchronization occurs, messages are constructed and exchanged in each of 3 coordinate directions (with barrier synchronizations between each message exchange event), and after all messages have been received, a barrier preceding the next time step is reached.  This is not by any means an optimal implementation.  However, if the message passing were overlapped with computation, this template code would lose its potential usefulness to the DoE. as a benchmark of SMP cluster network performance as well as computational power.  This benchmark aspect of the template code proved very useful.  It highlighted bugs in the prerelease MPI message passing on the Power Challenge Array which stimulated improvements in those libraries.
This 3-D template code is set up to solve a “real” problem of considerable interest to our DoE. sponsors, namely, the nonlinear unstable growth of perturbations to the shape of an interface between two gases of different densities which is accelerated impulsively by the passage of a strong shock front.  To obtain performance numbers for the Power Challenge Array, we used the Power Challenge machine from ARL in our lab together with 3 Power Onyx machines loaned by SGI, all interconnected by a HiPPI network with one HiPPI channel per machine.  Each machine had 2 GByte memory, and each Power Onyx had 8 R-8000 CPUs running at 90 MHz.  To balance this computational power, we needed to use 10 of the 75 MHz CPUs of the machine provided by ARL in order to simulate a machine with 8 processors running at 90 MHz.  We believe that the performance figures listed below should be very close to those we would actually obtain on an array of 4 identical 8-processor machines.  In the table below we treat the data as if it came from such a configuration, but the reader should be aware that this is not precisely the configuration that was actually used.  Also, in this table the speed-ups and efficiencies are based upon a comparison with an impossibly fast computation in which each cell in each 1-D sweep is updated at the speed of the sPPM 1-D module code.  In a full 3-D computation, y- and z-sweeps are less efficient, because of the need to transpose the fundamental data used by the scheme.  The cost of data copying made necessary by these transpositions and multitasking is not included in these estimates either.  Our data concerning the speed differentials for the individual 1-D sweeps indicates that these computational costs which have not been included in the reference single-processor timing would increase that timing by at least 8%.
Performance of the 3-D sPPM  ASCI  Template Code
on a Simulated SGI Power Challenge Array with 4x8 MIPS R-8000 CPU’s @ 90 MHz
(Efficiency and Speed-Up Estimated from Isolated 1-D Module Speed:  Underestimate)
GRID	Computation Time (sec)	Message Passing Time (sec)	Speed-Up (Underestimate)	Efficiency (Underestimate)	Gflop/s	
128*128*128	3.98	2.05	15.2	47.6%	1.70	
256x256x256	29.07	2.66	23.1	72.3%	2.59	
512x512x512	219.22	6.61	26.0	81.3%	2.91	

IIId.   The 3-D sPPM Template Code with Asynchronous Communication
The above performance data convinced us that despite our ability to send extremely large messages over the HiPPI channel at up to 85 MB/s using our own implementation of MPI, until something is done to reduce the cost of processor barrier synchronizations and explicit interactions with the operating system associated with message passing, we must drive the cost of message passing down by overlapping this function with computation.  This is much more difficult to do in 3-D than in 2-D, unless the grid in one of the three dimensions is very coarse, as it is in climate and ocean modeling.  The strategy which we are employing is still based upon the 8x8xN pencil of grid cells as the fundamental unit of work.  As in our tiled 2-D template code, we update these grid pencils first around the edges of the subdomain on our particular SMP.  This allows us to dispatch the messages from 4 of the subdomain’s 6 faces very early on in the calculation.  This procedure gives these messages a very long time to reach their destinations before the data they contain is needed to enable further computation.  The grid pencils in the interior of our subvolume contain data at each end of the updated pencil which must be sent to neighboring SMPs.  This data is dispatched as soon as the pencil has been updated.  Before beginning each 1-D sweep, we insert a barrier synchronization, which simplifies the logical flow of the program.  At this barrier we must wait for the last messages sent to arrive, but these messages are relatively small and cannot take very much time away from the calculation.
High computational efficiency can be obtained with this strategy so long as the message passing latency is of the order of 0.1 sec and the message passing bandwidth is of the order of 5 MB/s.  For this reason, we have chosen this template code as our vehicle for experimentation transcontinental distributed computing.  The performance of this template code on the ARL machine in the LCSE is given below.  Thanks to an SGI equipment loan, this machine now has 8 R-8000 CPUs running at 90 MHz and 10 running at 75 MHz.  The self-scheduled nature of the shared memory multitasking in this code makes possible the very high efficiencies recorded below for this heterogeneous collection of processors.  Clearly, this approach to parallel programming will yield benefits in a time sharing context as well, where differences in processor availability play the role of differences in processor clock speed in this test.  We also note that 82% efficiency on 16 processors for a small problem requiring only 5.5 seconds per time step is quite acceptable parallel performance for practical production work.
Performance of the 3-D PPM  Template Code with Asynchronous Message Passing
on a 128x128x64 grid
on the SGI Power Challenge Machine with 18 MIPS R-8000 CPU’s
( 8 @ 90 MHz and 10 @ 75 MHz )
#  CPU’s	Speed-Up	Efficiency	Sec / (2  t)	
1	1	100%	134	
2	1.97	98%	67.9	
4	3.91	98%	34.2	
8	7.69	96%	17.4	
16	12.0	82%	11.1	
18	11.7	72%	11.4	

Performance of the 3-D PPM  Template Code with Asynchronous Message Passing
on a 256x256x128 grid
on the SGI Power Challenge Machine with 18 MIPS R-8000 CPU’s
( 8 @ 90 MHz and 10 @ 75 MHz )
#  CPU’s	Speed-Up	Efficiency	Sec / (2  t)	
1	1	100%	1010.0	
2	1.98	99%	510.4	
4	3.94	98%	256.2	
8	7.84	98%	128.9	
16	14.4	98%	70.3	
18	15.8	97%	63.8	

IIIe.   The Red-Black Multigrid Relaxation Template Code
As part of the course which I teach in computational methods, I chose to develop an efficient multiprocessed implementation of the red-black relaxation scheme to solve the classic Dirichlet problem.  Essentially the same scheme can be used to solve the more general Poisson problem.  My original intent was to have the students parallelize this relaxation scheme for the SGI Power Onyx.  However, this proved rather involved, so I ended up doing this myself, writing up a detailed description of the implementation I developed and why I wrote it the way I did for the SGI machine.  I asked the students to explore the multiprocessing performance of this code and to report their findings and conclusions.  As a result, I have an additional template code, with documentation, which addresses a rather different class of applications than the sPPM templates discussed above.  The overall performance of this red-black relaxation scheme on a single 90 MHz R-8000 CPU is about 120 Mflop/s, while the compiler output indicates that the kernel loop should run at 167 Mflop/s.  I have compared notes with a collaborator at the University of Chicago who has worked with these multigrid relaxation schemes on the IBM SP-2.  He remarked that this code may be useful beyond its original pedagogical purpose.  We plan to work on this code together in order to get it into a more instructive template form which can be used in research applications and which can be parallelized across an entire Power Challenge Array rather than just on a single machine, as was done for my course.  This joint work will be part of a new Grand Challenge team project on stellar magnetoconvection which has recently been funded by NASA.
IIIf.   The  PPMLIB  Project
The PPMLIB project involves the development for widespread use of a package of library routines for computational fluid dynamics.  The library, PPMLIB, is specifically targeted for use on parallel computing platforms.  In the original DoE proposal nearly three years ago, two such platforms were specified, the Cray T3D and the Thinking Machines CM-5 Connection Machine.  In the early months of this grant, Thinking Machines declared bankruptcy, making the proposed library for the CM-5 a product of very limited usefulness.  After consultation with DoE., the strategy for producing PPMLIB was therefore altered.  A strategy was developed to target clusters of shared memory multiprocessors (SMPs), and to consider the Cray T3D/T3E and IBM SP-2 as degenerate cases in which the SMPs in the cluster each have only a single processor.  This strategy has involved extra work, but it was felt that this work would pay off in the generality of the resulting code.
During the last year, this strategy has proven to be precisely suited to the dominant emerging direction in the architecture of high performance computers.  The ARL Power Challenge Array project is built upon this vision of the high performance computing future, as is the program at NCSA.  So also are the ASCI “Blue” projects at Los Alamos and at Livermore.  These DoE labs have actually written algorithms from the PPMLIB project as “mandatory requirements” into their $100,000,000 ASCI Blue equipment procurement, which is now in process.  The Advanced Computing Lab at Los Alamos has also fallen in with this line, adding further legitimacy to our view of the computer architecture future.  As a final significant development, the recent Silicon Graphics merger with Cray Research makes virtually inevitable the ultimate transition of the entire high performance computing industry to the SMP cluster / Distributed Shared Memory paradigm.  It therefore would seem that our PPMLIB project is on the right track.
Central to the strategy of the PPMLIB project is the decision to support only a single Fortran‑P source version of the library routines.  Multiple target architectures and machines will be supported by multiple back-ends to our Fortran-P precompiler.  This strategy avoids the difficulty of supporting many individual versions of each library routine.  The manpower allocated to this project does not permit such a conventional strategy of individual code versions to be followed.  Because we must have only one source code, it is essential to determine the properties of that code, that is, the definition of the Fortran-P language extensions and coding style, before proceeding to develop the library.  This work has been largely completed over the last year.  We have worked extensively with the small SMP cluster of Silicon Graphics equipment in the LCSE, combining the ARL machine with other equipment loaned to us by SGI.  Our collaboration with SGI on the PowerWall visualization project brought to the LCSE during the fall months 2 additional 8-processor Power Onyx machines, each with 2 GB memory.  The NSF grant supporting our PowerWall project has recently purchased an Onyx with 4 R-10000 processors, and SGI is currently loaning to us 8 R-8000 processors, which have allowed us to augment the ARL machine to 18 processors.  Since all these machines are dedicated to development and experimentation rather than to production, the LCSE has been a unique environment for the software development which is critical to the success of this project.  Nowhere to my knowledge, except at the SGI home office in Mountain View, can one find an equally powerful SMP cluster configuration for algorithm and system software development and beta testing.
Over the last year, our algorithm testing has focused on the 4 code templates discussed above.  I refer to these as templates, because they represent different restructurings of our hydrodynamical codes, while the details of how individual grid cells are updated have no impact on these restructuring schemes.  In fact, our work with SMP clusters over the last year has led us to the conclusion that the difference schemes with which individual grid cells are updated may be bundled as Fortran-77 libraries appropriate for use on individual microprocessors, while all the shared memory multitasking and message passing between SMPs in a cluster can be handled separately in template codes.  The template codes which we have developed over the last year have been constructed from low-level multitasking primitives, which are equivalent to functionality found within the Posix threads standard, and from a small subset of the MPI message passing library.  We will rewrite these templates in Fortran-P, which should make them much more readable and easily understood by others.  We have made our templates available to collaborators around the country over the last year, before our Fortran-P precompiler is available.  Our intent has been to assist the community in understanding the principal programming techniques for SMP clusters in advance of our precompiler’s appearance.  Parallel programming experts can thus obtain the benefit of our early and extensive experience with this new architecture immediately.  Once our Fortran-P precompiler is ready, less anxious or less expert users will have the means to express highly efficient code procedures in as convenient and familiar a form as we are able to generate.  The theme, as from the start, is to invent no new syntax or semantics which is not absolutely necessary and to keep a strict focus on self-similar algorithms in order to make the task achievable in a simple and clean programming paradigm.
From the above discussion, the main outline of our strategy emerges.  We will bundle into PPMLIB two distinct sets of entities.  The first is a set, described in the original DoE proposal, of high quality numerical algorithms which operate on small strips, patches, or bricks of a logically regular mesh.  This first set is familiar;  it is a library of highly optimized Fortran-77 routines to which the global structure of the computation, or of the grid upon which it is carried out, is irrelevant.  These routines must be invoked in parallel by some other program structure in order to produce a parallel computation (although most of them are highly vectorized).  The complexities of the parallel machine architecture are shielded from these library routines, except insofar as they are optimized for cache memory microprocessors.  These routines can be constructed in very nearly optimal form in Fortran-77 independent of the particular processor upon which they will be executed.  In some cases, we will include in PPMLIB routines which are vectorized for Cray C-90, J-90, and T-90 type processors.  However, not all routines will have such counterparts, since these kinds of processors no longer dominate the high performance computing arena.  At this writing, the superscalar microprocessor chips which favor vector code, like the MIPS R-8000 and the IBM RIOS-2, are rapidly disappearing from the scene, to be replaced by chips, like the MIPS R‑10000, the DEC Alpha, and the PowerPC, which are designed with scalar code in mind.  Therefore, we are focusing our efforts principally upon this exciting new chip generation, which yields near vector performance on well structured scalar code.
The small library components, written in Fortran-77, are intended as atomic work units, that is, as uninterruptible tasks for single processors, in self-scheduled multitasked code.  If invoked properly, they will fit entirely into a normal-sized cache memory and execute extremely efficiently on most any microprocessor produced today.  These modules will be documented, and information about the work space they require will be made manifest, so that intelligent programmers can assure that these work atoms fit into the cache available on a given computing platform.  They will also, as indicated in the original DoE proposal, be bundled into larger, less flexible units in which these adjustments will have been made with care for the particular target machine.  The Cray T3E and the IBM SP-2 remain, as promised to the DoE sponsors, targets for PPMLIB, although the DoE’s ASCI Blue SMP clusters will be targeted as well (with our detailed work carried out on SGI machines, including the one from ARL, in our lab).
The larger PPMLIB bundlings, which can almost be thought of as complete codes, will follow the format of a set of distinct code templates, which we now intend to make available in Fortran-P source code format.  These code templates are far more general than the specific bundlings for hydrodynamics which we intend to offer to the community.  By making them available in Fortran-P source code, users should be able easily to adapt them to other purposes, much as we intend that they will adapt the smaller routines in PPMLIB to new uses by invoking them in new combinations and contexts.  It is our intent to offer a small number of choices for both 2-D and 3-D code templates.  Those which we have developed and tested to date were described earlier.
We are also planning to bundle as library modules single SMP tasks for an implementation on an SMP cluster.  If one builds a program for a single SMP, these tasks, which update large chunks of the regular grid, would play the role of whole programs.  They will contain within them, transparent to the larger program which invokes these library modules, self-scheduled loops over the atomic work units bundled as single CPU tasks in the lower-level PPMLIB routines.  The new feature which these single SMP tasks will introduce is task-to-task message passing via the intermediary mechanism of a shared memory (which need not be implemented in hardware).  As with the single CPU tasks, these tasks for SMPs are not explicitly aware of the existence of other tasks, but they nevertheless communicate information to those tasks by writing data into shared memory.  For the SMP tasks, the implementation of this process is straightforward for a distributed shared memory (DSM) machine.  For a cluster of SMPs, the implementation will involve MPI message passing.  This message passing, in the cluster case, will require knowledge, passed to the routine through its arguments, of the particular SMP which owns the portion of shared memory involved.  These SMP task elements of PPMLIB will therefore be the building blocks for asynchronous, self-scheduled, shared memory multitasking over the SMP cluster, whether the hardware supports such a global shared memory concept or not.  It is a key feature of this approach that the use of a global shared memory as an intermediate “message board” for the SMP tasks makes it possible to assign any of these tasks to any SMP on the network without having to inform that SMP of the identity of any other SMP which is executing any other task.  This essential feature, which comes at the cost of doubling both the memory and network bandwidth requirements for the job, allows a very straightforward and flexible dynamical load balancing technique to be employed in the application code.
IV.   Transcontinental Distributed Computing
Last spring, our group at the LCSE decided, with the encouragement of members of the DoE’s ASCI project and Larry Smarr of NCSA, to embark upon an experiment in transcontinental distributed computing.  We had realized that such a feat was possible for years, if only network bandwidth comparable to FDDI could be provided.  After all, our original experiment at SGI in Mountain View on the billion-cell turbulence calculation had used FDDI networking and involved message passing events separated by minutes (so that latencies of seconds could be tolerated).  This concept of transcontinental distributed computing became more interesting once the NSF set up its initial vBNS network, offering clear ATM OC-3 channels to facilities it interconnected.  The DoE plans a similar, and perhaps superior, network to join its 3 defense programs labs for the purposes of its ASCI project.  Through the NSF’s vBNS connections program, we expect to tap into both the NSF and DoE fast research networks at a sustained bandwidth of 10 MB/s from the LCSE.  NSF is supplementing one of our grants so that such a connection serving the LCSE along with other meritorious research projects at the University of Minnesota can be part of NSF’s beta test program for the vBNS.  Thus we hope to have this network connectivity in operation by October or November of 1996.
Our group at the LCSE took the lead in this transcontinental distributed computing project in the spring of 1995 in entering the HPC Challenge for the Supercomputing ’95 conference in San Diego.  Our collaborators consisted of Michael Brown and Brian Totty at SGI in Mountain View, John Towns at NCSA, and Charles Nietubicz at ARL.  We felt that the largest part of the work on this project was to get an appropriate code running.  The first 2 template codes discussed above were chosen, and the 3-D sPPM template with asynchronous message passing was developed (but not completed) with this project expressly in mind.  Other team members were repeatedly invited to develop other suitable application codes, but none were ultimately offered to us.
In November, our SGI collaborators in Mountain View arranged for our team to obtain exclusive access to a system of three 6x16 Power Challenge Arrays, to be interconnected with FDDI rings to approximate transcontinental computing conditions on the so-called “I-WAY.”  We hoped that this system, completely under our control, could give us the opportunity for a dry run.  Over the Thanksgiving weekend, Kevin Edgar from the LCSE went out to SGI in Mountain View to work on this dry run.  He was supported by Michael Brown and Brian Totty in California and by Steven Anderson and Joe Habermann in Minnesota.  This dry run helped us to identify a number of important bugs in SGI’s still unreleased implementation of MPI, but it did not give us the performance benchmarks which we desired or the scientific results which we had planned for.  After this dry run in Mountain View, Steve Anderson and Joe Habermann kept the project going with help mainly from John Towns at NCSA, who obtained the necessary permissions for us to use the I-WAY (these had inexplicably not been granted).  Although much progress was made and much of value was learned, this project proved to be too much of a stretch for the HPC Challenge in 1995.  Nevertheless, there is intense interest in the community for this project and, indeed, for the concept of transcontinental distributed computing for tightly coupled applications which it promises.

Performance of the Interstate Distributed Computing Experiment.
7 MIPS R-8000 CPUs at the LCSE in Minnesota
and 7 more at NCSA in Illinois, Connected via the Internet.

In the continued absence of the planned OC-3 link to Chicago, we decided to exploit the DS-3 link which already exists and which delivers bandwidths ranging from 200 to 600 KByte/sec during the night for the right kinds of messages sent between the LCSE and NCSA in Illinois.  We used the tiled 2-D sPPM template code described earlier to perform a single gas dynamics problem using an equal number of processors in each location.  Two grids were used, namely, 768x256 and 1536x512.  The conditions of this experiment were far from ideal, but the results, which are shown in the figure above, are encouraging nevertheless.
The 2-D grids used in the above experiment in distributed computing were rather small, and in order to achieve high latency tolerance, many subtiles, and hence significant redundant computation, must be involved.  We should expect better results, even for this simple Internet connection, for the directionally split 2-D PPM code discussed earlier.  For this code, the need for redundant computation in order to tolerate network latency and bandwidth limitations is almost entirely removed.  Again, two grid sizes were used, 1024x1024 and 2048x2048.  For the larger problem, using seven 90-MHz CPU’s in the LCSE and seven more in a machine at NCSA, performance was only 9% worse than that obtained by connecting two such machines with HiPPI in the same room.  These results are tabulated below.  For the larger problem, the superlinear speed-ups are obviously due to cache memory effects.

Performance of the Directionally Split 2-D PPM Template Code
on MIPS R-8000 CPU’s @ 90 MHz

GRID	#  CPUs	CPUs / SMP	Network Link (Bandwidth)	Speed-Up	Efficiency	
1024 x 1024	1	1		1	100%	
	2	2	1.2  GB/s	1.85	92%	
	4	4	1.2  GB/s	3.44	86%	
	7	7	1.2  GB/s	6.16	88%	
	14	7	100  MB/s	11.94	85%	
	14	7	Internet	6.69	48%	
2048 x 2048	1	1		1	100%	
	2	2	1.2  GB/s	2.26	113%	
	4	4	1.2  GB/s	4.38	109%	
	7	7	1.2  GB/s	7.56	108%	
	14	7	100  MB/s	14.9	106%	
	14	7	Internet	13.0	93%	

It is useful to estimate the network bandwidth which would be required to support this same sort of transcontinental distributed computing with the 3-D PPM code.  An sPPM computation on a 256x256x256 grid, which is carried out at 1 Gflop/s and 92% efficiency on a single 12-processor SMP, can be performed at each node of a distributed network in order to compute a larger problem.  Consider two such SMPs working together to update a grid of 512x256x256 cells.  The PPM computation is performed in 3 separate directional sweeps.  In 2 of these, all necessary messages can be dispatched very early in the computation.  In the remaining sweep, the grid is updated one pencil at a time, where each pencil has 8x8x256 cells on each network node.  Each pencil update takes 0.182 sec, after which a message of 4.61 KByte must be dispatched to the other SMP.  We require enough network bandwidth per processor in the SMP for this message to be transmitted in the 0.182 sec interval during which the next pencil is updated.  In this case, the required bandwidth is  12x0.00461/0.172 MByte/sec = 0.322 MByte/sec.  Note that in this estimate a 10 msec (interstate) message passing latency has been included.  Since there are 85.33 bundles of 12 pencils to be processed in this directional sweep, a sustained network bandwidth of 322 KByte/sec will make this computation very efficient indeed, despite tremendous geographical separation of the computing resources involved.  In practice, our group prefers to use the more accurate PPM scheme, which takes 60% more computing time to update each cell and which must afterward pass messages which are 7/3 longer.  For this scheme, a sustained network bandwidth of 470 KByte/sec is required to process the 512x256x256 grid cooperatively with two 12-processor SMPs running at 1 Gflop/s each.  We plan to test these projections using Silicon Graphics SMPs at the LCSE and at NCSA as soon as they are interconnected over the vBNS.
The example of a 512x256x256 PPM gas dynamics computation just discussed illustrates the principles involved in distributed computing of tightly coupled problems.  However, the aggregate computing resource assembled is only 2 Gflop/s.  This much computing power is already available at the LCSE in two Silicon Graphics SMPs.  Far greater computing power is presently available at ARL and at NCSA.  If, using such resources, 4 times the computational power which was assumed in the above example were applied to a single problem, roughly 2 MByte/sec of sustained network bandwidth would be required for such cooperative computation.  In this case, an aggregate grid of 512x512x512 could be updated in the same time as required for the 512x256x256 grid of the previous example.
There is another important factor to consider here.  An important reason to combine computing resources over a network is to get exploratory jobs done very quickly, perhaps even interactively.  As we reduce the problem size by a factor of 2 in each dimension, each pencil update takes only half as long, but the message which must be sent over the network during this time remains the same size.  Thus, for a computation on a 256x256x256 grid, which would take 8 times less time per mesh update as the previous example of a 512x512x512 grid, a sustained network bandwidth of 8 MByte/sec would be required.  A vBNS connection between the LCSE at the University of Minnesota and ARL or NCSA at ATM OC-3 bandwidth would meet this requirement, even if only half of this bandwidth were available for high performance networking research experiments.  Through the planned vBNS connection between the LCSE and NCSA, we will have the opportunity to demonstrate this exciting capability.
It should be noted that the problem sizes discussed here are not unreasonably large.  The largest PPM computation to date was performed in September of 1993 on an experimental Silicon Graphics computing system consisting of 16 Challenge XL machines with 20 processors each.  These machines were interconnected on 20 FDDI rings, and the PPM code implementation used in that experiment on a 1024x1024x1024 grid could tolerate message passing latencies of seconds.

V.   Preliminary Evaluation of the MIPS R-10000 Processor Performance
In early April, with the assistance of SGI, the LCSE was able to obtain the first 4-processor R-10000 board sold in the state of Minnesota.  This board was installed in the Onyx machine which drives our PowerWall display.  Because of the early shipment of this board, the second-level caches of the processors were only 1 MByte in size rather than the customary 4 MBytes, and we were told that other minor features of these processor chips were not as they were intended ultimately to be in quantity commercial shipments.  We also used the standard R-8000 compiler in our testing with these processors.  Thus we regard our performance results as preliminary.  Our R-10000 processors were clocked at 190 MHz and had peak theoretical performance of 380 Mflop/s.  During this same interval we had access, through Los Alamos National Laboratory, to 300 MHz DEC Alpha processors, with 4 MByte off-chip caches, the same chips which, we were informed, would be used in the Cray T3E (but, of course, without the off-chip caches).  These chips, along with the new HP chip, to which we did not have access, constitute a new wave of high performance processors optimized for the execution of scalar code.  If indeed these new chips can match the R-8000 performance on vector code, at a roughly two times lower cost per processor, and also deliver nearly the same performance on scalar code, there may no longer be a reason to use inherently inefficient numerical algorithms just because they vectorize.  This could usher in a new age of development of significantly more powerful numerical algorithms.  Our performance tests, reported below, indicate that this new age may well be at hand.
We carried out several tests of the R-10000, but the most relevant set, and one which duplicates most of the important results of our other tests, is our study of 4 different PPM or PPM-like gas dynamics codes.  These codes do not execute the same numerical algorithms, but all the algorithms are members of a family with many common characteristics.  The most important differences for this study were in the Fortran-77 implementations of these algorithms.  These include two highly vectorized implementations, the sPPM and PPM14 code kernels, and two highly scalar implementations, the PPM19 and MUSCL16 code kernels.  We used the Cray C-90 as a flop counter and as a means of measuring the degree of vectorization of the source code.  We gave the Cray the identical code versions which are optimized for microprocessors with caches, so the sPPM and PPM14 codes did not execute at their true potential on the Cray, which is about 450 Mflop/s for both codes.  However, our study was not about the C-90, since the characteristics of this processor are well known to all.  We also note that the C-90 compiler attempted to vectorize inherently scalar loops in the PPM19 code.  These loops, presumably, were vectorized by means of scatter/gather primitives, and, naturally, performance was poor.  The point here is to note that the vectorizability of this PPM19 code kernel, reported by the Cray compiler to be 93.3%, is clearly an overestimate in the normal sense of the term “vectorizability.”  The truth of the matter is better represented by the Mflop/s reported by the Cray’s hardware performance monitor for this code.  Most of the loops in PPM19 contain if-tests which, under certain conditions, cause control to jump to the end of the loop.  PPM19 is coded aggressively, so that a minimum of unnecessary flops which are normally performed in vector implementations (such as PPM14) are performed here.  PPM19 is an inherently more complicated algorithm than PPM14, and it yields more accurate results;  however, it performs only 1715 flops to update a grid cell in 2-D, while PPM14 performs 2671.  The 50% excess flop count is therefore a slight underestimate of the price in unnecessary computation which is paid for vectorization of even a totally regular gas dynamics problem.  This cost in flops should be kept firmly in mind when interpreting the results in the tables below.

Characteristics of Four Gas Dynamics Code Kernels
as Measured by the Cray C-90 Hardware Performance Monitor

CODE	CPU	Flops / Cell	% VECTOR	seconds per time step	Mflop/s	
sPPM	C-90	1647	99.8%	1.19	363	
PPM14	C-90	2671	86.3%	4.28	164	
PPM19	C-90	1715	93.3%	3.91	115	
MUSCL16	C-90	950	0.1%	15.7	15.9	

In the table above, the entries which are related directly to “time-to-solution” for this hydrodynamics problem on a 512x512 grid are those in the column labeled “seconds per time step.”  For the codes PPM14 and PPM19, which are fairly close in an algorithmic sense (although PPM19 is better) but fairly far apart in implementation, these entries are directly proportional to the time to update the grid, and hence to the overall cost of a computation on that grid.  The entries for sPPM and MUSCL16 in this column cannot be directly compared in this way, because these are both very simplified algorithms, and they do not deliver the same quality of result on a given grid as the PPM codes, or, for that matter, as each other.  The way in which MUSCL16 was coded is worth explaining.  It was derived directly from a completely vectorized code.  This code was transformed from a series of vector loops over a single strip of the grid to a single scalar loop over this strip which performs identical arithmetic but which makes no reference to temporary vector variables.  This transformation, which can be performed for any of the PPM family of codes, was done over 2 years ago with the hope of squeezing the data context for the grid cell update into the tiny on-chip cache memory of the DEC Alpha processor used in the Cray T3D.  The expected performance speed-up on that chip was not observed, but the code transformation is interesting because it results in what is likely to be the near absolute minimum cache memory requirement for the algorithm, at the cost of a total scalarization of this vector algorithm.  The vector parallelism of the algorithm (remember that the algorithm is left absolutely intact, identical to its vector implementation) is totally masked from the vectorizing compiler.  Before this set of performance tests was performed, one might have doubted that any code could bring the C-90 performance down to the 15.9 Mflop/s which appears in the table above.
Now that the characteristics of these gas dynamics code kernels have been explained, the performance figures for the MIPS R-8000, MIPS R-10000, and DEC Alpha processors can be interpreted.

Characteristics of Four Gas Dynamics Code Kernels
as Measured on the 90 MHz MIPS R-8000 (4 MB cache),
the 190 MHz MIPS R-10000 (1 MB cache),
and the 300 MHz DEC Alpha (4 MB cache)
updating a 512x512 grid

CODE	CPU	Peak Mflop/s	seconds per time step (Real*8)	seconds per time step (Real*4)	Mflop/s
(single
precision)	
sPPM	R-8000	360	6.87	3.89	111	
PPM14			7.32	5.89	119	
PPM19			8.47	7.45	60.3	
MUSCL16			4.87	3.84	64.9	
sPPM	R-10000	380	6.00	4.58	94.3	
PPM14			6.65	4.87	144	
PPM19			5.68	4.03	112	
MUSCL16			3.68	2.72	91.6	
sPPM	DEC Alpha	600	7.10	5.77	74.8	
PPM14			9.20	7.13	98.2	
PPM19			6.52	5.03	89.4	
MUSCL16			4.81	3.96	62.9	

One may draw several interesting conclusions from this table.  First, it is clear that the R-8000 chip, with its software-pipelined multiple floating point units, prefers vector code.  However, only a factor of 2 better performance is obtained with vector code on this chip, unlike the C-90, which delivers over a factor of 20 improvement in going from scalar to vector implementations.  Nevertheless, time to solution is minimized for the vector version of PPM on the R-8000.  On the R-10000 and DEC Alpha chips, with only a single set of floating point units per processor, there is very little penalty for scalar code, and, consequently, time to solution is minimized for the scalar PPM.  On these chips, there is also no discernible penalty, at least with present compilers, for minimizing the memory footprint of the algorithm by completely giving up any possible vectorization benefits.  This could mean that the DEC Alpha chip, implemented in the Cray T3E with no off-chip cache, could maintain the performance tabulated above, at least when the scalarizing code transformation used in MUSCL16 is invoked.  The good scalar performance of the new processor chips means that a wide variety of essentially unvectorizable numerical algorithms which perform very few gratuitous flops could now burst onto the scene, making these computers even more productive than previously expected.
As indicated in the table, the R-10000 is delivering roughly one useful floating point result on every second clock tick on scalar code.  It would be difficult to improve very much on this performance by building true vector capabilities into such a processor.  Cray XMP performance for our tuned implementations yields about one useful flop per clock tick.  For this vector processor, memory bandwidth is fully matched to the processor capabilities.  The R-10000 has only a third of the XMP memory bandwidth, measured in words per clock, and this bandwidth is delivered only from its cache, not from the main memory.  Modifying this processor chip to deliver full pipelined performance of one flop per clock on vector code might not be cost effective, when compared with simply coordinating the efforts of two unmodified R‑10000s sharing a common memory.
We also note that the tabulated speeds of the DEC Alpha chip are a much smaller fraction of this processor’s theoretical peak.  The reason for this is a mystery.  This table thus should serve to caution prospective computer buyers to examine more than the peak performance of the products.
One final comment should be made.  The brave new world of efficient scalar programs is not, necessarily, a world of fast scalar programs.  For programs to run fast, they must be both efficient on a single CPU and be parallel, with high parallel efficiency.  From the many results reported here, it is clear that vectorizable, regular computations can run with high parallel efficiency, both within a single SMP and across several SMPs.  For scalar programs, we expect high parallel efficiency to be achievable, but with greater difficulty.  The new difficulty which these programs will pose is one of dynamic load balancing.  Within a single SMP, the self-scheduled, barrier-free, multitasking programming style which we have employed in our template codes and which will be facilitated by our Fortran-P precompiler will solve this problem.  Over a cluster of SMPs, where this shared memory multitasking cannot be used in the most straightforward way due to insufficient network bandwidth, dynamic load balancing will have to become a new and critically important part of the program.  This function, which will involve dynamic reassignment of grid subdomains from one SMP on the network to another, with some mechanism to avoid potential thrashing, is closely akin to the normal function of an operating system.  As the discussion of our template codes, in this report and in the attached documents, should make evident, fast execution of a program on a single processor already requires a programmer to set up his or her own scratch variable arrays, which are resident in cache, with explicit code to copy data in and out of these arrays from and to main memory.  Thus programmers must already implement caches in their codes if they insist on performance.  The dynamic load balancing sections of future parallel programs will, in a similar fashion, be a mirror of operating system functions, effectively guiding these systems in performing their tasks effectively.
It is a tantalizing prospect to consider actually receiving assistance in dynamic load balancing from a network operating system.  Just which hints to such a system would be necessary within the program remain to be determined.  Programmers will need to decide in any case how their jobs can be decomposed into small parts which can be moved from SMP to SMP but which are not so small that unacceptable overheads result in moving them.  This requirement is similar to the need for programmers to determine the proper size of individual processor tasks, or work atoms, for efficient multitasking within an SMP.  Perhaps an ability to invoke certain operating system functions from within a parallel program would ease the programming difficulty here.  It would be nice to simply identify the subdomains and the code which operates upon them and to leave the movement of these units of work over the network, with all the complexities of boundary message passing which that movement implies, to the operating system.  This great day might be a while in coming, and therefore it may prove worthwhile in a future Fortran-P precompiler to assist programmers in this dynamic load balancing task.  The restriction of Fortran-P to self-similar algorithms implies that all message passing has a tremendous regularity.  The treatment within a Fortran-P program of all cells as interior cells also means that subdomains could, in principle, be merged together in a straightforward fashion to avoid the overhead of redundant computation in overlapping fake boundary cells on a single SMP.  For self-similar algorithms, an effective dynamic load balancing program thus may have a universal character, suitable for support from a precompiler.  Clearly, this is a subject for much further research.
VI.   High Resolution Visualization Systems and Software
With major funding from the National Science Foundation and new collaborations with the National Center for Supercomputing Applications (NCSA) and the Electronic Visualization Laboratory (EVL), our group at the LCSE has extended its early work on high resolution visualization systems and software.  We have purchased a PowerWall display system consisting of 4 Electrohome 9500 video projectors, each with 1600x1200 pixel resolution.  During the coming year, we hope to build this display into a wall of our lab, so that its use will be more convenient.  The PowerWall was originally developed together with SGI in 1994 with support from ARO through the AHPCRC and from the DoE., NSF, NASA, and ONR through research grants to our group.  The first PowerWall display system was demonstrated in the SGI booth at Supercomputing ‘94 in Washington, D.C.  It was driven by 2 Power Onyx machines, each with 2 Reality Engines-2, and by a farm of 24 Ciprico disk arrays, delivering sustained data rates in excess of 300 MByte/sec.  With the help of exceptional discounts from SGI, we have purchased an Onyx machine to drive the PowerWall display with dual “Infinite Reality” engines, each with 2 Raster Managers.  Thus we are able to drive all 4 high resolution panels of the PowerWall from a single Onyx chassis.  This configuration effectively halves the principal system cost for the PowerWall.  Two years ago, we consulted those at SGI who were finalizing the Infinite Reality engine specifications and made sure that this more cost-effective system design for the PowerWall would be possible.  Our Onyx, an upgrade of the Power Onyx which SGI loaned to our lab during the last year, has 2 GBytes of main memory and 4 R‑10000 processors.  It also has the maximum possible texture memory, thanks to assistance from Silicon Graphics.  Its 256 MBytes of texture memory should enable truly interactive parallel volume rendering of data sets from some of our largest simulations.  Using new DoE funding we are planning to attach to this system 4 Fibre Channel RAID systems, with an aggregate capacity of 280 GBytes and sustained data transfer rates of 320 MByte/sec.  The 4 RAID controllers are being loaned to the LCSE by our industrial partner Ciprico.
Now that our PowerWall system in the LCSE is coming together, we are beginning the software development needed to turn this display system into a more generally useful visualization tool.  We have continued our collaboration with Silicon Graphics in this enterprise, and we have added collaborators at NCSA, at the EVL, and at JPL, where PowerWall systems have recently been installed.  Porting of CAVE applications to the PowerWall immediately expanded the software base for this new system.  We have also built a PowerWall version of our group’s popular volume visualization software, in various forms.  This work has involved collaborations with a group in Norway and with Mike Palmer at Caltech.  Mike also worked closely with Brian Totty at SGI in Mountain View over the last year in speeding up his parallel adaptation of our group’s perspective volume renderer.
At Minnesota, our team has built specialized software which enables the generation of smooth movies from a relatively simple and interactive specification of a viewing path.  Graduate student Wayne Peterson used this software to create our first such movie with the data being viewed evolving in time as the view position changed.  We showed this movie at Supercomputing ’95 in the SGI booth and on the PowerWall (renamed by the conference organizers as the “NII Wall”) on the floor above the exhibit hall.  Wayne used this work as the basis for his Masters project in our department of computer science.  The PowerWall has proven to be a surprisingly popular visualization system, considering its cost.  SGI has received several orders for PowerWall systems, and at their request we are working to find a mechanism for integrating these systems through an appropriate third party arrangement with our industrial partner GeneSys.  GeneSys plans to provide full support to commercial customers for PowerWall software developed in collaboration with the LCSE.  It will also provide PowerWall installation services.  We note that SGI has included a description of the PowerWall in its recent product brochure.  To reduce the cost of an “entry level” PowerWall system, we have been developing specialized software which allows a PowerWall display to be driven for movie playback by 4 low-cost SGI systems.  Through the LCSE’s support arrangement with SGI, we hope to obtain the loan of 4 of SGI’s new Indigo systems for experimentation with this low cost method of driving the PowerWall and for further related software development.
VII.   High Performance I/O and Storage Subsystems and Software
Through efforts led by Tom Ruwart and augmented by Matt O’Keefe and his team, our group in the LCSE continues to work with state-of-the-art disk and tape I/O systems for use in visualizing and analyzing extremely large data sets.  Traditionally, these data sets have come from our group’s large fluid flow simulations with our PPM gas dynamics codes.  These data sets continue to drive our I/O system developments, but to these applications we have recently added a few different major data-intensive applications in the LCSE.  These new applications deal with data which is in digital form but which has not been generated by computer simulation.  Instead, these data sets come from experimental observations.  The first of these is work in the group of George Wilcox, of the department of pharmacology.  This work attempts to understand the neurological structure of the brain.  Voluminous microscopy data of rat and mouse brains is studied with interactive volume rendering techniques developed in our group.  The second major data intensive application is the Automated Plate Scanner project of Roberta Humphreys in the astronomy department.  In this project the Palomar Observatory’s all-sky survey has been digitized using a very high resolution plate scanner.  In this application, one is not only interested in interactive visualization of the data but also in automated recognition of special objects and sorting of these objects into classification bins.  Through our present and planned NASA projects, we expect this management of 2-D data from observing stations to become a larger focus of our I/O and storage work in the LCSE.  We expect to benefit from specialized software for such applications which will be developed at JPL for use with their PowerWall system.  We are advising JPL on storage subsystem purchases and expect that they will duplicate our experimental systems as their system software matures.
VIIa.   Fibre Channel Disk Subsystems
One of the LCSE’s industrial sponsors is Seagate Technology, a local manufacturer of high performance disk drives with whom our research group has worked on collaborative projects for many years.  Seagate is donating to the LCSE 36 disk drives, each with 9 GB capacity, incorporating that company’s latest Fibre Channel interface technology.  These drives have dual 100 MByte/sec Fibre Channel interfaces right on the drives.  We are working with Seagate to take advantage of these powerful new interfaces to build more capable and more cost effective RAID systems.  Over the next year, we will be building disk array systems based upon the Fibre Channel Arbitrated Loop (FCAL) technology which should offer unprecedented performance at unprecedented price points.  With industrial sponsor Ciprico we hope to collaborate on schemes to stripe over these FCAL systems to build extremely fast and cost effective disk subsystems.  Our initial goals are driven by the 300 MByte/sec data appetite of the PowerWall display, but we see no reason to stop at this particular figure.  GByte/sec disk subsystems should be a possibility with SGI’s new generation of computer systems.  We expect that in about 18 months, SuperHiPPI networking will begin to become available, with GByte/sec bandwidths matched to the transfer rates of the disk subsystems which we intend to build together with our industrial partners.
We have already had early Fibre Channel disk drives from Seagate in the LCSE for testing and software development.  Our initial work involved disks which were demonstrated by Seagate Technology at COMDEX last fall.  Because no 100 MByte/sec Fibre Channel interface boards were then available for SGI systems, we used a Pentium PC from Gateway 2000 for this demonstration.  However, we had several visitors in the LCSE from Prisa, a small company based in San Diego which is building Fibre Channel boards for SGI systems, working with our equipment to test early interface boards.  We are now collaborating with Prisa and SGI, and Prisa is becoming a formal sponsor of the LCSE.  As a start, Prisa has given us very special pricing on its early Fibre Channel boards for use in our Onyx, thus effectively donating 2 of the 4 boards to our lab in support of this effort.
VIIb.   Network Attached Disk Subsystems
One of the prime advantages of the Fibre Channel disk systems we intend to build at the LCSE is their potential to act as network attached devices.  This means that although requests for data must be processed by a host CPU on the network, the actual data stream delivered to the client need not pass through that CPU.  This arrangement frees up expensive CPU resources on the network and avoids network congestion.  LCSE Assistant Director Tom Ruwart tested this concept several years ago in a collaboration with Control Data Systems.  The result of this collaboration, our group’s SGI 300-series machine acting as the host to a striped disk subsystem delivering 10 MByte/sec network data transfers, was demonstrated in the AHPCRC booth at Supercomputing ’92 in Minneapolis.
Silicon Graphics took this network attached disk idea and pursued it, with the result that SGI demonstrated in its booth at Supercomputing ’95 a half-terabyte disk subsystem serving data to our PowerWall display at about 100 MByte/sec over 4 HiPPI channels through a network switch.  We see this technology demonstration as a major achievement, which will affect data-intensive computing in many application areas.  We believe that the Fibre Channel technology which we are pursuing in the LCSE offers still further improvements in performance and reductions in cost.  For a site like ARL, Livermore, or Los Alamos, where both classified and unclassified computing need to be accommodated, the network attached disk concept should have a special usefulness.
It should also be noted here that we have worked with SGI and LCSE industrial sponsor MTI (Micro Technology Inc., formerly National Peripherals) to demonstrate a terabyte file system at Supercomputing ’95.  This file system was constructed from MTI’s RAID-5 systems attached to 2 SGI Onyx machines which were used to drive the NII Wall (a big PowerWall) at the conference.  (The NII Wall was a joint project of the LCSE with the EVL and NCSA.)  Actually, the formatted capacity of this disk subsystem turned out to be just under a terabyte.  However, we reconstructed this disk system in the LCSE for a technology demonstration in mid December, and we added LCSE disk arrays to push the formatted capacity over the terabyte boundary, which seems to our collaborators to have some special significance.  The significance to us is to demonstrate the enabling technology, both software and hardware, which makes possible such cost effective storage systems.  Disk reads at 520 MByte/sec and disk writes at 200 MByte/sec were performed from the SGI Power Challenge machine provided by ARL.  The disk system was therefore fast as well as large.  Encouraged by this success, MTI and GeneSys are working with us to develop fast Fibre Channel disk subsystems from Seagate’s new disk drives.
Our work with network attached disks in the LCSE, led by Tom Ruwart, will involve the Ph.D. thesis project of Steve Soltis, under the supervision of Matt O’Keefe in the Department of Electrical Engineering.  This project, which involves building mechanisms and supporting system software to use a network attached array of disks to implement a sort of extended shared memory system, is outlined in Steve Soltis’s project proposal.  This proposal formed the basis of Steve’s oral exam, which he passed this spring.
VIIc.   High Performance and High Capacity Tape Subsystems
Over the last few years, Tom Ruwart has been leading a collaboration with Ampex to integrate their high-speed, high-capacity D-2 tape drives with computing systems from Silicon Graphics.  This work has also involved Ciprico in SCSI driver development.  A result of this work, which was begun 4 years ago in the AHPCRC, has been cost effective data archiving accessible to individual departments or projects which must deal with unusually high data volumes.  Last summer, our group acted as the catalyst for the coordinated acquisition of Ampex model 310 D-2 tape units by the LCSE, members of our turbulence grand challenge team at the University of Colorado, Boulder, and at Lockheed, and by 2 of the NSF supercomputer centers which we use, the Pittsburgh Supercomputer Center (PSC) and the Cornell Theory Center (CTC).  NCSA has expressed interest in joining this consortium of Ampex D-2 users by upgrading its present Ampex D-2 unit to the new model.  ARL may also wish to watch this project closely, since it has already invested in the older D-2 technology.
In addition to the integration of the Ampex D-2 drive with the SGI machines, our group, working with GeneSys, has generated a layer of software which creates a virtual tape drive.  As with any helical scan tape technology, it is of prime importance to keep the data streaming to the tape.  Interruptions in the writing process caused, for example, by the overhead of packing many small files into a single  tar  file on the tape can cause significant performance degradation.  Such interruptions of the data stream also place extra stress upon the mechanical tape transport system.  We have therefore generated an intermediate software layer which overcomes these difficulties by writing the data first to a “virtual tape drive” and only subsequently streaming this data to the physical tape.  For our applications, the 15 MByte/sec transfer rates of these Ampex tape drives are an essential reason for their purchase, and hence our development of this software.  Last winter, on the coldest day in Minneapolis history, we demonstrated this software to a group of visitors from Ampex in California.  They were sufficiently impressed to loan to the LCSE a model 410 D-2 system, which includes a robotized stack of up to 7 tape cartridges, with over a terabyte of on-line storage.  We are using this system for software development, and are making it available to companies in the area, such as Control Data Systems, for testing.  This system allows multiple users to be served without operator intervention, which is of great benefit to us.  We are developing software for analysis and visualization of very large data sets which relies on this system to replace a much more expensive large disk farm.  For sequential operations, such as movie generation, this has proven to be very effective and convenient, freeing up precious high speed disk resources to interactive tasks.  We are discussing with Ampex the possible loan to the LCSE of a model 810, with multiple tape drives and very large on-line storage capacity.  Industrial partner GeneSys, our collaborator in these software development projects, may become involved in supporting these tape systems for Ampex customers.
VIId.   I/O and Storage Technology Demonstration at the LCSE
Our team’s work with I/O and storage systems was demonstrated in the LCSE on December 18-21, 1995.  Working with our industrial partners, Silicon Graphics, MTI, and GeneSys, we demonstrated the terabyte file system, our integration of Ampex and Silicon Graphics systems, and the PowerWall to a number of industrial visitors as well as the general public.  The systems demonstrated were shown at the Supercomputing ’95 conference.  Nevertheless, interested companies, in some cases, flew in representatives for this more private technology viewing and discussion.  This technology open house was described in the St. Paul paper on the front page of the business section, with a color picture of the PowerWall.  This very positive report prompted over 100 phone calls for information and gave the University of Minnesota some much appreciated good publicity.  We plan to repeat this event next year, but not so close to the Christmas holiday.
VIII.   Summary
This first year of LCSE operation has been very productive.  Activities in our group funded by the NSF, the DoE, NASA, and ONR have grown in intensity, in size, and in productivity as a result of this new facility, new grants, and the direct involvement of our industrial partners SGI, MTI, Ciprico, Seagate Technology, GeneSys, Prisa, DEC networking, and Ampex.  Our work on I/O and storage systems as well as on high resolution visualization systems and software proceeds at a vigorous pace.  We are making periodic technology demonstrations, such as those at the Supercomputing ’95 conference in San Diego, which can serve to introduce these developments to the supercomputing community.  Finally, we are pursuing a project in transcontinental distributed computing for relatively tightly coupled applications.  We expect to be connected to the NSF’s high speed research network, the vBNS, and will use this connectivity to demonstrate this new capability to the community.  This network connection could provide a means of carrying out distributed computing experiments with NCSA, as well as with Los Alamos, Livermore, and possibly ARL.  The same software which makes these demonstrations possible will also permit us to use the very fast disk subsystems which we are developing at the LCSE together with the SGI Power Challenge machine provided to us through our joint research project with ARL and SGI to perform unprecedented out-of-core computations for the numerical simulation of turbulent fluid flows on very fine grids.