A recently introduced general interpolation method based on reproducing kernel Hilbert space (RKHS) theory has been quite successful in constructing a number of potential energy surfaces. A straightforward implementation is slow when large numbers of ab initio points are involved, since the computation time is proportional to the number of points. The algorithm introduced here precomputes and stores the sums necessary to generate the surface, allowing the surface computation time to be independent of the number of points. The method is tested on its ability to reproduce a prior global analytic fit to the potential surface for the N(D-2)+H-2 system based on a 7X7X7 grid of points. The RKHS interpolated surface is found to exactly reproduce the 343 points on which it is based, and has a root mean square (rms) error of 14.2 kJ mol(-1) elsewhere, while the prior analytic fit has a rms error of 25.1 kJ mol(-1) at the points used to fit the surface. With a 16X16X16 grid as input the RKHS surface had a rms error of 1.32 kJ mol(-1). The interpolated surface with the new algorithm can also be evaluated twice as fast as the prior analytic fit.