Barycentric interpolation is a numerically stable algorithm to evaluate Lagrange interpolating polynomials. An implementation in Julia is available at BarycentricInterpolation.jl.

The procedure takes as input a set of nodes X_k = \{ x_1, x_2, \ldots x_k\} and initially computes a set of weights W_k = \{ w_1, w_2, \ldots w_k \}.

Then, evaluating the k polynomials at the point x involves:

(This is pseudocode adopted from the above mentioned package just for description purposes)

```
for i = 1:k
vals[i] = Wn[i] / (x - Xn[i])
end
```

(The complete algorithm checks for the case `x == Xn[i]`

and then sets `vals[i] = 1.0`

and all other entries of `vals`

to `0.0`

. Also, there is a final normalization step.)

I was wondering if there is a “Barycentric algorithm” that deals with the case when `x::Interval`

(from IntervalArithmetic.jl).

Specifically, if `x`

contains `Xn[i]`

then the output bounds become arbitrarily loose `[0,+/-∞]`

. Any inputs in dealing with this will be very helpful!

(It’s a little amusingly ironic that the barycentric algorithm is used because it is numerically robust and avoids cancellation errors, yet seems to not work immediately with interval arithmetic.)