Walking that back a bit. Given a set of points (x[n],y[n]), the interpolated value can be calculated in O(N^2) time if the x coordinates change, and in O(N) tim

Walking that back a bit. Given a set of points (x[n],y[n]), the interpolated value can be calculated in O(N^2) time if the x coordinates change, and in O(N) time if only the y values change. For common cases, this means O(N^2) initialization, and O(N) evaluation, even if the Y values change. Source: https://people.maths.ox.ac.uk/trefethen/barycentric.pdf https://people.maths.ox.ac.uk/trefethen/barycentric.pdf This paper is (admittedly) a significant improvement to the method I've been using to date. Probably the most significant improvement is a method for calculating Lagrange polynomials that provide stability guarantees. After having read it, I'm going to point you to that paper, instead of describing how I did it. The application that I'm familiar with, is interpolation of audio signals, and fractional delay lines. In these cases, x[N] can be considered to have the fixed values [0,1,2,3,...], and then subsequent interpolations can be calculated using an offset into the buffer of audio samples. Initialization takes O(N^2). Evaluation takes O(N).