Computational science research and development

The mission of CSS is to help realize the end-to-end scientific simulation environment envisioned by the NCAR Strategic Plan. The specific role of CSS is to develop much of the critical software and intellectual infrastructure needed to achieve the plan's ambitious goals. To more accurately reflect this mission, the Computational Science Section has been reorganized in 2002 into three functional groups: Mathematical Methods, High Performance Technology Tracking, and the Earth System Modeling Framework Project. This alignment enables CSS to more effectively track computer technology, learn to extract performance from it, pioneer new and efficient numerical methods, create software frameworks to facilitate scientific advancement -- particularly interdisciplinary geoscience collaborations -- and share the resulting software and findings with the community through publications, talks and websites.

CSS has been both productive in research activities and successful in obtaining grant money in 2002. The section published six papers, including one review article, in peer-reviewed journals in 2002 and has six more currently in press. CSS had five successful proposals in 2002, bringing into the section a total of approximately $3,500,000 of new funding, and $160,000 of new equipment.

The largest single funded activity -- for $2,800,000 over three years -- is the Earth System Modeling Framework project. The ESMF project is a highly visible interagency effort to build model creation and model-coupling software infrastructure for the Earth science modeling community. ESMF is well underway: in the past year it has successfully met the first four NASA milestones of the project, and is currently in the design phase of the software system.

Another significant CSS initiative in 2002 was the CSS-led effort that resulted in NCAR becoming a member of the Gelato Federation, an open software consortium focused on the Itanium architecture. This resulted in CSS and NETS Web100 funding, in the form of both salary for software development and equipment that amounted to approximately $360,000. This funding will both accelerate development of the Spectral Toolkit and provide NCAR with visibility within the Linux community.

As a result of this progress, CSS has also grown dramatically over the past year. Eight new staff members have been added:

• Within ESMF, four full-time software developers, one software integrator, and a half-time administrative assistant have fully staffed ESMF as the NASA-funded project gets underway.
• The scientific staff in CSS has grown from one to four in 2002. The CSS core mathematics and computer science research program has been expanded by the reclassification of one software engineer IV (Dr. Steve Thomas) to Scientist III, the addition of a Scientist I-level mathematician, an expert in conservative advection, and a joint appointment at the Scientist II level with the Computer Science Department at the University of Colorado.
• The mathematical libraries effort in CSS has added another software engineer to augment the Spectral Toolkit Mathematical Software Library under the HP-sponsored Gelato Federation Initiative.

High performance technology tracking

CSS technology-tracking efforts in 2002 have centered on three activities:

• Understanding and evaluating the new IA-64 (Itanium) computer architecture
• Converting the 36-processor Compaq AlphaServer named prospect into a functioning Linux cluster renamed lhotse
• Performing benchmark tests and developing performance models in support of the ARCS procurement

The Itanium's very long instruction word (VLIW) architecture represents a dramatic departure from the traditional superscalar RISC microprocessor and CISC-like Pentium architectures used in the geosciences departments at most universities today. As Itanium microprocessors become plentiful in the geosciences community, access to open-source high-performance mathematical libraries optimized for this architecture will become important to scientific progress. It is clear that, to be well positioned for future advances, CSS needs to learn to performance-program the Itanium processor. To this end, in the fall of 2001, CSS acquired a four-processor Rx4610 Itanium-I server from Hewlett Packard. This system was evaluated through synthetic benchmarks and also by porting several NCAR applications, including CCM3.

This activity also included an inter-comparison of the usability and robustness of the Linux and HP-UX operating systems on the introductory Itanium system and extensive evaluation of early releases of the Intel and HP C++ and Fortran-90 compilers for Itanium. This work was carefully documented and presented to HP engineers in Fort Collins, Colorado in February 2002. A few months later, HP agreed to sponsor NCAR as a member of the Gelato Federation, an organization devoted to the advancement of the Linux/Itanium environment and the distribution of open source software. As a member, CSS will receive $200,000 over two years to fund a software engineer to help optimize the STK libraries and other applications for the Itanium architecture. In May 2002, CSS also responded to and won an HP RFP for free Itanium-2 hardware: six Itanium II systems, two dual-processor servers, and four workstations worth $130,000 were awarded. Two of these will be provided to NETS to conduct Web100 development and testing in collaboration with the Pittsburgh Supercomputing Center.

The second HPC technology-tracking activity has been the conversion of the largely unsuccessful $750,000 AlphaServer system prospect into a Linux system. Early in this effort, CSS established a collaborative relationship with a new Scientist I -- Peter Thorton in CGD -- who needed approximately four months of dedicated computing to perform a data-intensive, embarrassingly parallel calculation to run a 2,000-year spin-up calculation of a detailed CONUS land surface model at an unprecedented 1-km resolution. This application mapped ideally to the new Alpha-processor-based Linux cluster now renamed lhotse. This calculation, the first of its kind, has been successfully completed.

Mathematical methods research

The mathematical research effort within CSS is divided into two principal activities:

• Developing performance-portable, highly efficient, open-source numerical libraries for use by the mathematical, geosciences, and engineering communities
• Investigating mathematical techniques and software implementations for producing accurate, efficient and scalable numerical dynamical cores for general circulation and turbulence models

Numerical Library Development

The first version of the spectral toolkit (STK) was released in 2002. Initially, access was restricted to early adopters in the NCAR scientific community; currently the initial release of STK software and documentation may be downloaded from the CSS website http://www.scd.ucar.edu/css/software/stk/ The release version of STK provides support for complex and real one-dimensional FFTs. Although the software is written in C/C++, it supports a simple and easy-to-use Fortran-callable interface. All STK functions are threadable, and some, such as the two- and three-dimensional complex and real FFTs, have built-in parallelism implemented using pthreads. These parallel functions in the STK library have been applied and thoroughly tested in the construction of a parallel turbulence model in collaboration with NCAR scientists Jack Herring and Yoshi Kimura. This 3D model was run on parallel computers at resolutions up to 1024 points in each spatial dimension. A distributed-memory version of the multi-dimensional 3D FFT has achieved over 200 GFLOPS sustained performance on the Pittsburgh Terascale computing facility.

Paul Swarztrauber has developed new software that implements the computational methods for computing Gauss points and weights and for computing harmonic projections. The new method for computing Gauss points and weights has already been incorporated into the SHEREPACK subroutine gaqdi. This software will ultimately be integrated into the Spectral Toolkit.

The projection algorithms based on the weighted orthogonal complement method have also been incorporated into SPHEREPACK as subroutines shpei, shpe, shpgi, and shpg, and they perform projections on equally spaced and Gauss-distributed grids. The new SPHEREPACK library can be accessed and downloaded at http://www.scd.ucar.edu/css/software/spherepack/

Elliptic solvers for semi-implicit nonhydrostatic NWP models

The elliptic problems in semi-implicit nonhydrostatic atmospheric models are typically poorly conditioned, nonseparable, contain cross-derivative terms, and often are nonsymmetric. In 2002, the collaborative effort between NCAR researchers Steven Thomas (SCD), Joshua Hacker (MMM), Piotr Smolarkiewicz (MMM), and Roland Stull (University of British Columbia) has led to a class of effective Krylov methods, i.e. a conjugate residual (GCR) algorithm preconditioned with a 3D direct solver, using standard tridiagonal inversion in the vertical. They have developed a horizontal spectral preconditioner as an alternative to the more standard and much simpler line-Jacobi relaxation scheme. However, the efficacy of the spectral preconditioner requires neglecting the cross derivative terms and homogenization (e.g., averaging) metric coefficients over the computational domain. Because such a compromise causes a substantial departure of the preconditioner from the original elliptic operator, it is not obvious a priori whether it leads to a competitive solver. They evaluated the robustness of the proposed approach over a broad range of representative meteorological applications, and documented its superior performance in the context of a three-time-level semi-implicit semi-Lagrangian all-scale weather-prediction/research model.

Spectral Element Atmosphere Model (SEAM)

Further progress was made in FY2002 on both 2-D shallow water and 3-D primitive equations SEAM models based on a spectral element horizontal discretization. Support for hybrid coordinates was added to the primitive equations dynamical core in anticipation of introducing physics into the model. Studies were done on new filtering techniques, based on the work of Fischer and Tufo, that do not require communications.

The domain decomposition software in the SEAM spectral element model was improved in FY2002 in two important ways. To minimize communications costs of the pure MPI paradigm on SMP clusters, the package now uses a hierarchical domain decomposition strategy based on the METIS package from the University of Minnesota. The communications library first decomposes the element mesh on the cube sphere into node-domains based on the number of SMP nodes; then each node-domain is fed back into METIS to be decomposed into processor-domains.

The second improvement is the introduction of decompositions based on the nesting of so-called Peano or Hilbert space-filling curves. Such curves map the two-dimensional element grid on the cube sphere onto a one-dimension sequence of elements that is generated in a self-similar hierarchical way. Hilbert curves permit one to produce efficient load-balanced domains with quite good locality with a trivial domain decomposition algorithm. The package is currently limited to a number of elements along a cube edge factorable into powers of two and three. The space-filling curve decomposition, as applied to the cube grid in SEAM, is the subject of a paper in preparation by John Dennis in CSS.

The integration rate of the semi-implicit shallow-water model version of SEAM has been further improved by introducing an additive overlapping Schwarz preconditioner for the iterative conjugate gradient solver. The Schwarz preconditioner has a much smaller memory footprint compared to the block Jacobi preconditioner employed originally. The Schwarz preconditioner also has slightly better convergence properties than the block Jacobi. This is because it overlaps with neighboring elements, and so does a better job of communicating information about the solution throughout the system.

Finally, this year Stephen Thomas and Henry Tufo have been awarded $500,000 for three years from the NSF to explore adaptive mesh refinement strategies using the SEAM dynamical core. The funding includes support for a post-doctoral fellow, a graduate student research assistant, and annual workshops that will emphasize the dynamical core research within CSS. These workshops will offer financial support for young researchers and students to visit NCAR to gain experience working with high-performance numerical libraries, dynamical cores, test cases, and visualization tools.

Double Fourier dynamical core

Further progress on the double-Fourier approach included work on a semi-implicit semi-Lagrangian scheme with NCAR ASP post-doctoral fellow Anita Leighton. However, work on this aspect of the project has stopped within CSS due to the departure Dr. W. Spotz from the section. A presentation by Drs. Leighton and Spotz will be given at the 2002 PDEs on the Sphere meeting in Toronto, Canada.

In related work, Paul Swarztrauber worked with Bill Spotz to compare four candidate projection algorithms for use in the Double Fourier Method:

• The standard method using forward and backward Legendre transforms
• The direct method, which uses a single projection matrix when this approach results in fewer operations
• The fast multipole method
• The weighted orthogonal complement method

Swarztrauber and Spotz have shown that for spectral truncations up to N=200, the weighted orthogonal complement method has the lowest operation count, best cache utilization and best overall timings.

Cache-optimized spectral dynamical core (BOB)

Richard Loft and Leo Rivier in CSS have developed, in collaboration with Lorenzo Polvani at Columbia University, an efficient, cache-optimized spectral transform dynamical core for low-cost Linux Beowulf clusters. This code, called Built-on-Beowulf (BOB), outperforms CCM by substantial factors (3-5) at relevant climate resolutions on these low-cost systems.

The BOB code organizes a vorticity-divergence formulation of the shallow water equations (or primitive equations) along the lines of the method described by Rudi Jacob-Chien. This formulation permits one to write the Legendre transform as a blocked matrix-matrix multiply. The memory footprint of the Legendre transforms are made smaller by computing the associated Legendre polynomials "on the fly." This reduces memory traffic and increase cache reuse.

The Ph.D. thesis of Rivier was based on simulations of jet formation in the upper atmosphere of Jupiter using the BOB code. An article recently published in Mon. Wea. Rev. documented the design and performance of BOB. Subsequent to publication of this paper, there have been requests for access to the code from researchers at GFDL, MIT, and Texas A&M University.

Coronal Mass Ejection

Dr. B.C. Low of NCAR has proposed a theory for the formation of Coronal Mass Ejections (CMEs -- magnetic "bubbles") based on the forced (gravitation and pressure) and force-free elliptic equation for the solar magnetic stream function. The theory relies on the total magnetic energy of the associated magnetic field.

Dr. Steve Thomas of CSS, in collaboration with Dr. Natasha Flyer and Bengt Fomberg of the University of Colorado, have applied a variety of advanced numerical methods to solve this highly nonlinear elliptic problem. These include both second-order finite difference methods and a Fourier-Chebyshev pseudo-spectral method on the semi-infinite half-plane in spherical coordinates.

Because the potential far-field solution decays in inverse proportion to radius, Dr. Flyer's original contribution was to use an exponential conformal map for the problem. The far-field solution is then matched to the inner solution at an artificial outer boundary using a radiation boundary condition (RBC), derived from the eigen-functions of the linear far-field operator. The solutions of the elliptic problem depend on one critical parameter, and the team is currently exploring the use of arc-length continuation combined with Newton's iterative method of approximation to compute solutions as well as looking for bifurcation points.

The sought-after magnetic "bubble" solutions for the forced equation have been found, and initial results in the force-free case are very encouraging. The team believes that a nonlinear bifurcation may lead to the formation of a bubble, and we are tracking the magnetic field energy as a function of the parameter.

Earth System Modeling Framework

Computational Science Section 2003 activities

CSS will exploit the dramatic successes of FY2002 in FY2003 by:

• Meeting the ESMF Project milestones
• Evaluating and developing mathematical software for the Itanium architecture as Gelato Federation members
• Further developing the SEAM and BOB dynamical cores and related partnerships
• Starting the newly funded initiatives in the areas of ocean forecasting and spectral element adaptive mesh refinement research and development activities

During FY2003, the ESMF team anticipates completing the following activities which are required for the next four project milestones

• Installing a high-performance cluster at MIT
• Releasing a prototype API and software
• Achieving partial ESMF compliance for a subset of the ESMF testbed applications, which include the CCSM, WRF, the NCEP forecast suite, and other codes
• Making initial delivery of three interoperability experiments with new coupled configurations for testbed codes

The focus of the work with lhotse in FY2003 will center on installing hardware and software to create a scalable I/O subsystem for lhotse and installing the kernel modifications needed to support high-speed message passing using the Quadrics interconnect. This will allow lhotse to be used as a compute platform for the GTP program.

In FY2003, the Spectral Toolkit will be further expanded in several areas. Early in FY2003, CSS intends to release additional STK functionality including 2D and 3D multithreaded real and complex FFTs and spherical harmonic Legendre transforms on the Gaussian grid. The performance of the Legendre transform will be based on hand-optimized real-matrix complex-matrix multiply algorithms that have already demonstrated high fractions of peak on a variety of microprocessors. An effort to optimize all STK functions for the Itanium-2 processor, funded by Hewlett-Packard Corporation, will augment this effort.

Support for spherical harmonic Legendre transforms on equally spaced grids will be added early in 2003, as will additional operators, such as those for computing the vorticity or divergence from given input velocity fields and vice versa. The ultimate goal is to incorporate the essential functionality of SPHEREPACK into STK, including the new functions recently added by Paul Swarztrauber. Work will also begin next year developing 2D and 3D multi-threaded Spectral Element quadrature operators and the associated mathematical framework that will allow users to easily construct efficient multi-threaded and distributed-memory applications based on the spectral element method.

STK is envisioned as gradually forming the basic mathematical infrastructure for all of the dynamical core research in the section. This will provide a useful way to disseminate the core knowledge obtained by this research to a broad community of users that may in turn employ these components in ways not anticipated by the creators of the STK library.

The first step in this direction will be to develop a new implementation of the BOB dynamical core, which will employ the spherical harmonics and FFT functions in STK. Both 2D shallow-water and 3D primitive equations versions of BOB will be implemented on top of the STK libraries. We expect this core to surpass the performance of the original BOB implementation.

The 2D semi-implicit scheme employing the additive overlapping Schwarz preconditioner for spectral elements will be extended to the 3D primitive equations. An eigen-mode decomposition in the vertical permits one to solve a 2-D Helmholtz problem for the generalized pressure in each vertical layer. To further improve the integration rate acceleration factor over the explicit Eulerian formulation, we plan to introduce conservative semi-Lagrangian advection, based on the work of Dr. Nair, for the dynamics and constituent transport in the spectral element model.

To date, only the Held-Suarez idealized physical forcings and the baroclinic instability tests of Jablonski and Williamson have been performed on the spectral element core. The next step with the section's spectral element model will be to systematically introduce and test physics, or more precisely, introduce the core into the CAM2 framework of CCSM. Decompositions of the cube sphere element grid which load-balance the physics will have to be developed and tested. The incorporation of physics will likely follow the development and evaluation of several intermediate test cases, such as an adiabatic atmosphere with smooth orography, advection of passive tracers over orography, followed by test cases with active water vapor, such as the aqua-planet of Williamson.

Two dynamical cores and libraries workshops are planned over the next three years to transfer CSS research in these areas to the broader NCAR-university community. These workshops will be tutorial in nature with an emphasis on the construction of dynamical cores using CSS libraries. Test cases and visualization will also be emphasized, including the NCAR Graphics language NCL. Finally, CSS staff will continue to spread knowledge of their findings through presentations and publications in peer-reviewed journals.

In addition to these continuing activities, CSS plans to begin building a computer science program. To this end, in FY2003, CSS intends to pursue two proposals. The first is a midrange NSF ITR proposal to investigate the suitability of IBM Blue Gene/L and Blue Gene/C advanced MPP architectures to geoscience simulation. This proposal is being developed in collaboration with computer scientists at the University of Delaware, Lawrence Livermore National Laboratory, IBM T.G. Watson Research Center, the University of Colorado, and the Computational Science Section within SCD.

The second proposal will be submitted to NASA and will involve the development of a Grid testbed between the NCAR Mass Storage System (MSS) and a Pentium-4 Beowulf cluster at CU to produce a fully production-capable, Grid-enabled, bio-geo-chemistry modeling system. This initiative will build the resume of CSS in the area of Grid computing, which could then provide the basis for additional Grid proposals in the future.

Educational outreach

Cecelia DeLuca in CSS organized and led a hands-on tutorial session February 20, 2002 to introduce Earth System Modeling Framework (ESMF) participants to the framework development environment. The training session was held at the NCAR Corporate Training Center.

Steve Thomas co-organized, in collaboration with Annick Pouquet of the Geophysical Turbulence Program, a scientific conference entitled Adaptive and High-order Methods with Applications in Turbulence. This conference was held at NCAR in March 2002. The workshop lasted three days and had 15 invited speakers, with an additional four graduate student and post-doctoral presentations. The conference brought together approximately 30 scientists interested in the numerical simulation aspects of turbulence.

Community service activities

Cecelia DeLuca in CSS served as one of three technical managers of the Earth System Modeling Framework (ESMF), which is developing the software infrastructure for building and coupling climate models, weather models, and data-assimilation systems. Collaborators include SCD, CGD, and MMM at NCAR; NOAA GFDL and NCEP; MIT; the University of Michigan; DOE ANL and LANL, and NASA GSFC, DAO, and NSIPP. The ESMF is being developed jointly by the Earth system community as a resource that will benefit the community by increasing software reuse, interoperability, and performance portability. Cecelia organized the first ESMF community meeting, which was held in Washington, D.C. on May 30, 2002. The purpose of the meeting was to review a draft ESMF Requirements Document and elicit additional requirements from the broader community. More than 80 people attended the meeting. Additional input from the community was collected via online forms and archived on the ESMF website.

Cecelia DeLuca was reappointed CCSM Software Engineering Working Group (SEWG) Co-Chair for a second two-year term. She organized and participated in SEWG meetings at NCAR February 14-15, 2002 and at the CCSM workshop in Breckenridge, Colorado on June 25, 2002.

John Dennis presented a talk at the SEWG Breckenridge workshop on Trends in Computer Architectures. His presentation provided important technical background information to the attendees of the workshop.

Steve Thomas has refereed publications for Monthly Weather Review, Journal of Computational Physics, International Journal of Applied Mathematics, and Computer Science. Steve was elected to the Scientific Advisory Committee (SAC) for high-performance computing in Canada.

Richard Loft reviewed articles for the Journal of Computational Physics.

CSS was instrumental this year in bringing about UCAR membership in the Gelato Federation http://www.gelato.org/   The goals of the Gelato Federation are to provide open source software for the Linux/Itanium-2 computing platform. The architecture of the Itanium-2 represents a significant departure from previous microprocessor designs and holds significant promise as a key component in future high-performance computer designs. As a Gelato member, CSS has committed one software developer and has recently received $200,000 in matching funds from Hewlett-Packard, UCAR's Gelato sponsor, to fund a software engineer for at least two years. HP has also donated approximately $160,000 of hardware to support this effort. CSS's objectives as members of Gelato are three-fold: port and optimize CSS Spectral Toolkit mathematical libraries for Itanium-II, port and tune production NCAR Fortran codes to this environment, and improve, in collaboration with the NETS Web100 project, the TCP/IP network performance of these systems. These three activities will combine to enable UCAR researchers in the geosciences to better exploit the performance capabilities of this promising architecture.

Publications

Baker, A., J. Dennis, and E.R. Jessup, 2002: Algorithm Design for Fast Linear Solvers. Proceedings of VECPAR 2002, 5th International Meeting on High Performance Computing for Computational Science, Porto, Portugal, 429-442.

Dickinson, R.E., S.E. Zebiak, J.L. Anderson, M.L. Blackmon, C. DeLuca, T.F. Hogan*, T.M. Iredell*, M. Ji*, R.B. Rood*, M.J. Suarez*, and K.E. Taylor*, 2002: How Can We Advance Our Weather and Climate Models as a Community?, BAMS, 83.

Flyer, N. and P.N. Swarztrauber, 2002: The convergence of spectral and finite difference methods for initial-boundary value problems, SIAM J. Sci. Comput., 23, 1731-1751.

Loft, R.D., S.J. Thomas, J.M. Dennis, 2001: Terascale Spectral Element Dynamical Core for Atmospheric General Circulation Models. Proceedings of the Supercomputing 2001 Conference. ACM and IEEE Computer Society, CD-ROM.

Polvani, L.M., R.K. Scott, and S.J. Thomas, 2002: An initial-value problem for testing the dynamical cores of atmospheric general circulation models. Submitted to Mon. Wea. Rev.

Rivier, L., R.D. Loft, and L.M. Polvani, 2002: An efficient spectral dynamical core for distributed memory computers. Mon. Wea. Rev., 130, 1384-1396.

Spotz, W.F. and P.N. Swarztrauber, 2001: A performance comparison of associated Legendre projections. J. Comp. Phys., 168, 339-355.

Swarztrauber, P.N. and W.F. Spotz, Eds., 2002: Computing in Science and Engineering Special Issue on Climate Modeling. Comp. in Sci. and Eng.

Swarztrauber, P.N., 2002: On computing the points and weights for Gauss-Legendre quadrature. SIAM J. Sci. Comput., accepted for publication.

Thomas, S.J., R.D. Loft, and J. Dennis, 2002: Parallel Implementation Issues: Global vs. Local Methods. Submitted to Computation in Science and Engineering.

Thomas, S.J., J. Hacker, P.K. Smolarkiewicz, and R. Stull, 2002: Spectral preconditioners for nonhydrostatic atmospheric models. Submitted to Mon. Wea. Rev.

Thomas, S.J., R.D. Loft, 2002: Semi-Implicit Spectral Element Atmosphere Model, SIAM J. Sci. Comput., 17, 339-350.

Thomas, S.J., J.M. Dennis, H.M. Tufo, and P.F. Fischer: An overlapping Schwarz preconditioner for the cubed-sphere, submitted to SIAM J. Sci. Comput. special issue, Copper Mountain Conference on Iterative Methods, 2002.