Towards adaptive mesh PDE simulations on petascale computers

The advent of the age of petascale computing brings unprecedented opportunities for breakthroughs in scientific understanding and engineering innovation. However, the raw performance made available by petascale systems is by itself not sufficient to solve many challenging modeling and simulation problems. For example, the complexity of solving evolutionary partial differential equations scales as n at best, where n is the number of mesh points in each spatial direction.1 Thus, the three-orders-of-magnitude improvement in peak speed of supercomputers over the past dozen years has meant just a factor of 5.6 improvement in spatio-temporal resolution—not even three successive refinements of mesh size. For many problems of scientific and engineering interest, there is a desire to increase resolution of current simulations by several orders of magnitude. As just one example, current longrange climate models operate at O(300) km resolution; yet, O(10) km grid spacing is desirable to obtain predictions of surface temperature and precipitation in sufficient detail to analyze the regional and local implications of climate change, and to resolve oceanic mesoscale eddies [9]. Thus, although supercomputing performance has outpaced Moore’s Law over the past several decades due to increased concurrency [9], the curse of dimensionality imposes much slower scientific returns; e.g. improvement in mesh resolution grows at best with the onefourth power of peak performance in the case of evolutionary PDEs. The work requirements of scientific simulations typically scale as n. The power α can be reduced through the use of optimal solvers such as multigrid for PDEs and fast multipole for integral equations and N-body problems. Once α has been reduced as much as possible, further reductions in work can be achieved only by reducing n itself. This can be achieved in two ways: through the use of adaptive mesh refinement/coarsening (AMR) strategies, and by invoking higher order approximations (in space and time). The former place mesh points only where needed to resolve solution features, while the latter