A moving interpolation technique which provides an accurate representation of potential energy surfaces for polyatomic molecules is presented. The method uses the ab initio energy, energy gradient, and second derivatives calculated at dynamically important configurations. The interpolant of the energy and its derivatives converges to the exact value with increasing number of data. A procedure is given for finding the optimum configurations at which ab initio calculations are performed. The method is demonstrated by application to the six-dimensional surface of a diatomic plus diatomic reaction.