1. Overview
2. Linear Algebra for Dense Systems
   1. The BLAS
   2. LAPACK
   3. BLACS
   4. PBLAS
   5. ScaLAPACK
3. Sparse Linear Algebra
   1. The Sparse BLAS
   2. The PIM Library
4. Other Parallel Libraries
   1. PESSL
   2. NAG Parallel Libraries
   3. PETSc
5. Summary of Features
6. Parallel Libraries at the MHPCC
7. Example Programs
8. References and More Information

• Several parallel math libraries are available for the SP platform. Use of these libraries can reduce programmer development time, promote code portability and reuse, and improve code performance.

• Libraries may facilitate entry into parallel programming. Several libraries parallelize existing serial routines

  • Serial: example interfaces to LAPACK LU factor and solve routines:

     dgetrf(N, N, A, lda, pivots, ierror)
     dgetrs(trans, N, nrhs, A, lda, pivots, b, ldb, ierror)
     

  • Parallel: interfaces to corresponding ScaLAPACK routines:

     pdgetrf(N, N, A, ia, ja, desc_A, pivots, ierror)
     pdgetrs(trans,N,nrhs,A,ia,ja,desc_A,pivots,b,ib,jb,desc_b,ierror)
     

  Note: converting a serial program that calls LAPACK to a parallel program using ScaLAPACK requires additional work to manage the processors and data

• Current parallel math libraries will expand in functionality over time and others will become available

  • ScaLAPACK (Scalable LAPACK)

  • PESSL (IBM Parallel ESSL)

  • NAG Parallel Libraries

  • PIM (Parallel Iterative Methods)

  • PETSc (Portable, Extensible Toolkit for Scientific Computation)

• Some currently available functions are:

  • Linear and nonlinear matrix operations including parallelized BLAS

  • Direct and iterative linear solvers

  • Non-linear solvers

  • Eigensystem and singular value analyses

  • FFT (1- to 3-dimensional)

  • Random number generators

• This tutorial will present a brief overview of several libraries currently available.

  • Details of their use will not be given - it would require an individual tutorial for each library to accomplish this.

  • Examples and pointers to more information will be provided

  • The basis of LAPACK and ScaLAPACK will be presented to motivate the discussion

• BLAS (Basic Linear Algebra Subprograms) - set of Fortran 77 computational kernels for basic operations commonly used by linear algebra routines, de facto standards. 3 Levels of BLAS.

• Level 1 BLAS

  • In 1973, Hanson, Krogh and Lawson described the advantages of adopting a set of basic routines for problems in linear algebra [1]. Their specification was intended to:
    • Be the basic building blocks to aid in the design and development of numerical linear algebra software
    • Improve readability of code
    • Promote modular code that is more robust and maintainable
    • Promote efficiency (by use of best algorithms and by encouraging vendors to implement on their machines)
    • Improve portability - de facto standardization, machine dependencies can be "hidden" from higher level code

  • Primarily scalar and vector operations
    ex: z = ax + y, x and y are vectors and a is scalar

  • The Level 1 BLAS operate on an element at a time of a vector or matrix. The calculations are fine grained (the ratio of floating point operations to memory references is small).

  • Level 1 BLAS were very successful and have been used in a wide range of software including LINPACK.

  • See the Level 1 BLAS routine descriptions for details on the individual routines.

• Level 2 BLAS

  • Motivated by the development of vector-processing machines, Dongarra et al. proposed an extended set of BLAS in 1988 [2]. The level 2 BLAS were designed to perform linear algebra operations at a 'higher level' than the level 1 BLAS. Instead of simple vector operations, these routines are able to carry out matrix-vector operations, and to solve linear systems.

  • matrix-vector operations
    ex: z = aAx + by, x and y are vectors, a and b are scalars and A is a matrix

  • See the Level 2 BLAS routine descriptions for details on the individual routines.

• Level 3 BLAS

  • Proposed by Dongarra et al. in 1990. Targeted at matrix-matrix operations with the aim of providing more efficient, but portable, implementations of algorithms on high-performance computers, especially those with hierarchical memory and parallel processing capability.

  • matrix-matrix operations
    ex: D = aAB + bC, a and b are scalars, and A, B and C are matrices

  • See the Level 3 BLAS routine descriptions for details on the individual routines.

• BLAS Efficiencies
  • BLAS 1: Practical experience on emerging high performance vector, hierarchical memory and shared memory machines indicated that a coarser level of granularity was needed to obtain the higher performance offered by these machines.

  • The Level 2 and 3 BLAS improve efficiency on high performance machines by coarsening the granularity.

  • The level 2 and level 3 BLAS can be exploited by the use of block (or partitioned) algorithms

  • Block algorithms partition matrices into blocks and perform the work on each of the blocks

  • Block algorithms increase efficiency by
    • Control of data movement between levels of memory (cache/vector registers, main memory, disk memory)
    • Increase in reuse of data held in lowest level of memory (cache or vector registers)
    • Maximization of amount of data held in low-level memory

  BLAS Level	BLAS routine	O(Ops)	O(Mem Refs)	Ratio
  1	SAXPY z=ax+y	2n	3n	2:3
  2	SGEMV z=aAx+by	2n^2	n^2	2:1
  3	SGEMM D=aAB+bC	2n^3	4n^2	n:2

• LAPACK (Linear Algebra Package)

  • Successor to LINPACK and EISPACK incorporating the functionality of both

  • Portability addressed by use of BLAS as the building blocks

  • Efficiency gained by using block algorithms that exploit level 3 and level 2 BLAS

• LAPACK Routines - designed for dense or banded matrices, not for general sparse matrices

  • Systems of linear equations

  • Linear least squares problems

  • Eigenvalue and singular value problems, including generalized problems

  • Matrix factorizations: LU, Cholesky, QR, QZ, SVD, GSVD, Schur and generalized Schur

  • Condition and error estimates

• Target Machines

  • Vector processors

  • Shared-memory multi-processors

  • Workstations and PCs

• The LAPACK User's Guide is available from netlib at http://www.netlib.org/lapack/lug/lapack_lug.html

• BLACS (Basic Linear Algebra Communication Subroutines) - set of kernel communication routines for message passing parallel linear algebra routines with the added purposes to:

  • Aid design and coding

  • Promote efficiency by identifying frequently occurring operations of linear algebra
    • Operations are matrix based - readily transmit sub-matrices of a larger matrix
    • Include point-to-point communication, broadcasts, reduce operations and support routines
    • Context driven - safety mechanism to ensure that messages sent in one context cannot be received in another

  • Ensure portability via standardization of kernels

• Optimized version of the BLACS are available in several formats for various machines

• The BLACS are not a general purpose communications library (usually vector based rather than matrix)

• See the BLACS routine descriptions for details on the individual routines. Several BLACS program examples are also available.

• PBLAS (Parallel Basic Linear Algebra Subprograms)

• Functionality similar to BLAS - distributed forms of vector-vector, matrix-vector and matrix-matrix operations

• Simplify parallelization of dense linear algebra code - building blocks of parallel linear algebra. Allow the software writer to focus on performing message passing on subsections of matrices rather than at low level byte transfers.

• Designed so that local computations within a process are performed by the BLAS, while the communication operations are handled by the BLACS (actual implementations may/may not do this, however).

• Portability layer for ScaLAPACK - together with BLACS and BLAS - machine dependencies confined at this level, de facto standards

• See the PBLAS routine descriptions for details on the individual routines.

• A complete discussion (including examples) on PBLAS can be found at http://www.netlib.org/utk/papers/pblas/pblas.html.

• ScaLAPACK (Scalable LAPACK) is a library of high-performance linear algebra routines for distributed-memory message-passing MIMD computers and networks of workstations supporting PVM and/or MPI.

• Essentially a port of LAPACK - contains the same routines for solving systems of linear equations, least squares problems, and eigenvalue problems.

• Resulted from the ScaLAPACK Project involving Oak Ridge National Laboratory and several universities.

• Design goals:

  • Maintain performance (including scalability)

  • Retain portability

  • Stay as close as possible to LAPACK in calling sequence, storage, etc.

  • Promote modularity via set of linear algebra tools (BLAS, BLACS and PBLAS)

  • Use LAPACK algorithms when possible

  

• Performance issues

  • Data movement
    • Between vector registers/cache, memory and disk as for LAPACK

    • Between memory of distributed processors

    • Controlled via block algorithms that use Level 2 and 3 BLAS (as in LAPACK)

  • Load balancing to keep processors active achieved via 2-D block-cyclic data distribution
    • Load balance is improved with smaller block sizes

    • Communication costs are lower for larger block sizes since fewer messages of larger size are sent

    • There is a block size that balances the tradeoff between load balancing and message latency

• The ScaLAPACK Users' Guide is available from netlib at http://www.netlib.org/scalapack/slug/scalapack_slug.html

• The sparse BLAS are extensions to the BLAS

• Goals similar to those of BLAS:
  • Portability via de facto standardization

  • Efficiency via tuned implementations for various machines

  • Improved program clarity, modularity and maintainability

• Compressed storage of sparse vectors (non-zero entries)

• Sparse vector operations determined by examining existing sparse libraries for common operations
  • BLAS routines that perform non-vector operations (e.g. Givens rotations)

  • BLAS routines that operate correctly on compressed storage vectors (2-norm, sum of absolute values, constant times vector, etc.)

  • BLAS routines requiring sparse extensions

• Documentation, routines for a model implementation and more available from netlib at http://netlib2.cs.utk.edu/sparse-blas

• A set of sparse BLAS and a model implementation have recently been defined, no libraries based on them are yet available.

• PIM (Parallel Iterative Methods) is a library of FORTRAN 77 routines to solve systems of linear equations on parallel machines using iterative techniques.

• Main goals:
  • Allow user freedom to choose matrix storage, access and partitioning methods

  • Support portability across variety of parallel architectures and programming environments

• Mechanism to achieve goals: hide details of the following operations (user supplies the routines)
  • Matrix-vector product

  • Preconditioning step

  • Inner products and vector norms

• PIM components:
  • Conjugate-Gradients
  • Bi-Conjugate-Gradients
  • Conjugate-Gradients squared
  • Stabilized Bi-Conjugate-Gradients
  • Restarted stabilized Bi-Conjugate-Gradients
  • Restarted generalized minimal residual
  • Restarted generalized conjugate residual
  • Normal equation solvers
  • Quasi-minimal residual with coupled recurrences
  • Transpose-free quasi-minimal residual
  • Chebyshev acceleration

• For additional information, downloading software and the PIM User's Guide, see the PIM web page located at http://www.mat.ufrgs.br/pim-e.html

• PESSL (Parallel Engineering and Scientific Subroutine Library) is is a mathematical subroutine library from IBM designed to provide high performance for numerically intensive computing jobs running on IBM SP Systems or on clusters of IBM RS/6000 workstations. It is IBM's parallel analogue of its serial ESSL library.

• Includes subroutines for:

  • Linear algebra subprograms

  • Matrix operations

  • Linear algebraic equations

  • Eigensystems analysis

  • Fourier transforms (2 & 3D), convolutions and correlations, and related computations

  • Uniform random number generator

• PESSL components:
  • ESSL for core computational routines (over 400 subroutines)
    • BLAS - tuned for RS6000 processor
    • Suite of common numerical routines including a subset of LAPACK

  • BLACS - latest version uses the Parallel Environment (PE) Message Passing Interface (MPI) routines to communicate data between processor nodes and coordinate the execution of parallel tasks. Previous versions used MPL.

  • Subset of level 2 and 3 PBLAS

  • Subset of ScaLAPACK linear algebra, eigensystem and singular value analysis routines.

• Each subroutine offered in PESSL (version 1.2 +) has a related High Performance Fortran (HPF) subroutine. Provides the ability to make calls to PESSL from HPF applications.

• PostScript manual (over 700 pages!) available in /usr/lpp/pessl/doc/postscript on the MHPCC system.

• The IBM PESSL web page is located at: http://www.rs6000.ibm.com/software/sp_products/esslspec.html

• The NAG Parallel Libraries are the parallel analogue of the NAG Fortran Library. There is a PVM based library and an MPI based library.

• NAG library components:

  • NAG Fortran libraries for core computational routines

    • BLAS - the Fortran library contains untuned BLAS; tuned BLAS should be substituted when available

    • Suite of numerical routines including most of LAPACK

  • BLACS - communication is PVM or MPI based

  • Level 2 and 3 PBLAS

  • Subset of ScaLAPACK

  • Quadrature

  • Unconstrained minimization

  • Sparse linear algebra

  • Uniform random number generator

• PETSc (Portable, Extensible Toolkit for Scientific computation) is for both uni- and parallel-processor scientific computing:

• Especially intended for large-scale problems modeled by partial differential equations

• Contains suite of data structures and routines (linear and non-linear equations solvers)

• Uses MPI for message passing

• Provides interfaces to other libraries such as LAPACK and the BLAS

• Is available by anonymous ftp from info.mcs.anl.gov. Is supported, documented, and "will be actively enhanced for the next several years".

• PETSc components:
  • Vectors - set of level 1 BLAS-like serial and parallel vector routines

  • Matrices - set of routines and data structures for manipulating parallel sparse matrices

  • Krylov Solvers - set of Krylov space iterative solvers with or without pre-conditioners

  • Pre-conditioners - variety of pre-conditioners

  • SLES (Simplified Linear Equation Solver) - high level interface to linear equation solvers

  • SNES (Simplified Nonlinear Equation Solver) - routines to solve systems of nonlinear equations or unconstrained minimization problems

• The PETSc home page at Argonne National Lab is located at http://www.mcs.anl.gov/petsc

Areas of Numerical Analysis
Area	ScaLAPACK	PESSL	PETSc	PIM	NAG	AZTEC	SYISDA
Linear Equations Direct methods	yes	yes	serial	no	yes	no	no
Linear Equations Iterative methods	no	no	yes	yes	yes	yes	no
Least Squares	yes	yes	yes (iterative)	no	yes	no	no
Eigensystem Analysis	yes	yes	no	no	yes	no	yes
Singular Value Analysis	yes	yes	no	no	yes	no	no
Random Number Generation	no	yes	no	no	yes	no	no
Quadrature	no	no	no	no	yes	no	no
Fast Fourier Transforms (FFTs)	no	yes	no	no	no	no	no
Nonlinear Equations	no	no	yes	no	no	no	no
Unconstrained minimization	no	no	yes	no	yes	no	no
PDE	no	no	yes	no	no	no	no

Table of Basic Features
Feature	ScaLAPACK	PESSL	PETSc	PIM	NAG	AZTEC	SYISDA
Application Languages	F77	F77, C, C++	C, C++, F77	F77	F77	C, F77	C
Source Languages	F77	F77	C	F77	F77	C	C
Data types	Single precision, Double precision, Complex, Double complex	Double precision	Double precision, Double complex	Single precision, Double precision, Complex, Double complex	Single precision, Double precision	Double precision	Double precision
SP Communications Library	BLACS (MPI/MPL)	BLACS (MPI/MPL)	MPI	PVM, MPI	BLACS, PVM	MPI	PRISM, MPI
Availability	Public domain	Vendor	Public domain	Public domain	Vendor	Vendor with public research license	Public domain
Documentation	Technical reports & users guide	750+ pages	Draft manual	Manual included	Intro. and Tutorial	40 page manual	29 page User's Guide
Data structure	Traditional	Traditional	Transp. to user. Several formats available	Provided by user. Transp. to PIM	Traditional	Distrib. Modified Sparse Row, Distrib. Variable Block Row	?

 

Parallel Libraries at the MHPCC

General guidelines:

• Use LAPACK in favor of LINPACK

• Use tuned BLAS for the architecture you are running on when available. For SP/RS6K machines, link with the tuned ESSL library. The table below demonstrates the difference in performance between the BLAS library and ESSL which uses tuned BLAS.

Matrix Size	Routine	BLAS library		ESSL library
		Thin 64MB Node (Mflops)	Wide 256MB Node (Mflops)	Thin 64MB Node (Mflops)	Wide 256MB Node (Mflops)
500 order	SGEMM	158	154	224	224
	DGEMM	146	160	242	245
	CGEMM	50	51	203	219
	ZGEMM	64	64	219	241
1000 order	SGEMM	162	165	228	238
	DGEMM	147	160	166	227
	CGEMM	52	51	243	247
	ZGEMM	64	64	124	237

Note: Multiplies were of square matricies.C=C+A*B (except in ZGEMM because of space problems on thin nodes they were C=C+A*A)

• The PESSL library can be expected to give the best performance when it contains the routine you need and your problem size is optimal. The table below shows the effect of problem size on performance for the PESSL subroutine pdsyevx, which computes the eigenvalues/eigenvectors of a real, symmetric matrix. Note that the parallel speedup would be better for larger problems.

  problem-size 240x240			problem-size 480x480
  serial ESSL (sec)	#nodes	pESSL(sec)	serial ESSL(sec)	#nodes	pESSL(sec)
  1.5	1	2.4	11.2	1	9.8
  .	2	2.2	.	2	7.3
  .	4	1.6	.	4	4.9
  .	9	1.5	.	9	3.8
  .	16	1.5	.	16	3.5

Where they are:

• The public domain packages are under either /usr/nfs/packages/math (LAPACK, ScaLAPACK, PETSC) or /source/pd/math (PIM). Corresponding libraries are generally available in /usr/nfs/lib with associated include files in /usr/nfs/include. When available, man pages are placed in /usr/nfs/man.

• ESSL, PESSL: the libraries are in /usr/lib

• NAG libraries: /source/vendorcode/nag (except MPI Library)
  • Fortran Library (serial): /source/vendorcode/nag/nagfl16df

  • Parallel/PVM Library: /source/vendorcode/nag/nagfd01df

  • Parallel/MPI Library: /usr/nfs/packages/math/fmib201d

• LAPACK User's Guide
  (http://www.netlib.org/lapack/lug/lapack_lug.html)

• ScaLAPACK Home Page
  (http://www.netlib.org/scalapack/scalapack_home.html)

• PESSL: an on-line manual is available via InfoExplorer, set your DISPLAY variable then enter "info -l pessl" at the Unix prompt. A postscript version of this manual is available in /usr/lpp/pessl/doc/postscript. Caution: this manual is very large. You may want to select single pages for printing by using ghostview.

• NAG PVM library: a hypertext based on-line manual is available on tsunami.sp2.mhpcc.edu. To use the on-line manual, add the following lines to your .cshrc file:

       setenv EBTRC /source/vendorcode/nag/odunx01nd/data/.ebtrc
       alias dtext /source/vendorcode/nag/odunx01nd/bin/dtext
  

  Then enter dtext at the Unix prompt.

• PETSc User's Guide
  (http://www.mcs.anl.gov/petsc/manual.html)

• PIM: a postscript manual is available in /source/pd/math/pim20/doc

LU Factor and Solve

DESCRIPTION:

These programs demonstrate the use of the LU factor and solve routines to solve a dense system of linear equations. A serial example calling the routines from LAPACK is given together with a parallel example that using the ScaLAPACK library routines. The programs are written in Fortran 90.

FILES:

• Serial: ser_lu.f

• Parallel: par_lu.f

• Input: lu_solve.in

• Makefile

• LAPACK, ScaLAPACK, BLAS, PBLAS, BLACS
  NETLIB web server: (http://netlib2.cs.utk.edu/)

• PIM home page: (http://www.mat.ufrgs.br/pim-e.html)

• PETSc home page: http://www.mcs.anl.gov/petsc

• [1] Hanson R., Krogh F. and Lawson C. (1973) "A Proposal for Standard Linear Algebra Subprograms" ACM Signum Newsletter, 8, p. 16.

• [2] Dongarra J.J., DuCroz J., Hammarling S. and Hanson R.J. (1988) "An Extended Set of FORTRAN Basic Linear Algebra Subprograms" ACM Trans. Math. Softw., 14, pp. 1--32.

• [3] Dongarra J.J., DuCroz J., Duff I.S. and Hammarling S. (1990) "A Set of Level 3 Basic Linear Algebra Subprograms" ACM Trans. Math. Softw., 16, pp. 1--28.

• Notes from the "Workshop on the Use of Public Domain Software and Software Standards", in particular the contributions of Jack Dongarra and Sven Hammarling, conducted at the MHPCC.

• Pittsburgh Supercomputing Center's online software documentation for descriptions of the individual BLAS routines.

• Centre Interuniversitaire de Ressources Informatiques de Lorraine online software documentation for descriptions of the individual BLACS routines.

• Workshop material developed by the Albuquerque Resource Center.

• Freeman, T. L. and Phillips, C., "Parallel Numerical Algorithms", Prentice Hall, New York

• Golub, Gene and Ortega, James M., "Scientific Computing: An Introduction with Parallel Computing", Academic Press, Boston

• MHPCC tutorial first written by Lon Waters