The manipulation of two dimensional elevation maps is an important part of digital cartography. In many situations, these maps are computed by interpolating sparse data such as isolated elevation points obtained from stereo matching. In this paper, we present a surface interpolation algorithm based on variational splines which is well suited to massively parallel computers. Using multiresolution parallel relaxation, we can efficiently compute the interpolated surface and also have local control over its continuity and smoothness. We apply this technique to sparse elevation data and to elevation contours, and show how to add realistic fractal detail through stochastic relaxation. We also present a multiresolution decomposition algorithm and a fast parallel 3-D rendering algorithm.