barycentric-rational#

An interactive illustration of

M. S. Floater and K. Hormann, Barycentric rational interpolation with no poles and high rates of approximation, Numerische Mathematik 107 (2007) 315–331. doi:10.1007/s00211-007-0093-y

The method itself is a MATLAB-syntax script you can edit in the page. It runs in your browser through numbl; there is no server and nothing to install.

What the paper says#

Interpolating a function at a given set of points is easy to do badly. The degree-n polynomial through n+1 equally spaced points diverges as n grows, which is Runge's example. Classical rational interpolation, fitting a quotient p_M/q_N with M + N = n, often approximates better but offers no control over where the poles land, and they land inside the interval.

Floater and Hormann's construction is short. Fix an integer d with 0 ≤ d ≤ n. For each i let p_i be the polynomial of degree at most d through the d+1 points x_i, …, x_i+d, and blend those local polynomials together:

The results are that r has no poles anywhere on the real line for any d (Theorem 1), that its error is O(h^d+1) for d ≥ 1 whatever the node distribution, as long as f is smooth enough (Theorem 2), and that r can be rewritten in the barycentric form

with explicit weights (equation 18), which is cheap to evaluate. On a uniform mesh those weights are integers, nearly all equal, differing only near the two ends: 1, 4, 7, 8, …, 8, 7, 4, 1 for d = 3. That small change at the ends is what lifts the approximation order from O(h) to O(h⁴).

The d = 0 case is Berrut's earlier interpolant, and d = n is ordinary polynomial interpolation, so the family interpolates between the two.

What the page shows#

Four tabs, all driven by the same script:

• Interpolant — f, the rational interpolant r, and the nodes, with the degree-n polynomial and a clamped C² cubic spline as optional overlays, plus the pointwise error underneath.
• Blending & weights — the n−d+1 local polynomials and the normalised blending functions λ_i, which sum to 1 everywhere but have oscillating tails and no local support; and a stem plot of the barycentric weights, with the integers of Section 4 read off when the mesh is uniform.
• Poles — the denominator of r on the real line, drawn as a signed n-th root so that its sign and zeros survive the enormous dynamic range; all of its roots plotted in the complex plane, none of them touching the real axis; and the classical p_M/q_N alongside, with the poles it does put in the interval.
• Convergence — max error against n on log-log axes for a range of d, with the measured orders tabulated. With Runge's function on uniform nodes this reproduces Table 1 of the paper, and with the spline turned on, Tables 3 and 4.

The script#

The editor holds the whole method. The app only asks it for two functions, and optionally a third:

w = bary_weights(x, d)             % the weights, equation (18)
r = bary_eval(x, y, w, t)          % the barycentric form, equation (1)
[P, L] = local_blend(x, y, d, t)   % the blend of (4) and (5)   [optional]

Anything that satisfies that contract will drive all four tabs, which is the point of the alternative scripts in the method dropdown:

The last two are the counterpart to Schneider and Werner's theorem, quoted in the paper, that a pole-free barycentric rational interpolant must have weights that alternate in sign.

Running it locally#

npm install
npm run dev            # http://localhost:5173
npm run build          # type-check and bundle to dist/

Two test suites, neither of which needs a browser to be watched:

npm run test:matlab    # the .m layer against the paper's tables, via the numbl CLI
npm run test:matlab -- --full   # also n = 640
npm run test:browser   # the built app in headless Chrome

test:matlab runs the MATLAB layer outside the browser through a clone of numbl (set NUMBL_DIR; it defaults to ~/src/numbl) and checks it against the published numbers: the integer weight patterns of Section 4 for d = 0…4, the error columns of Tables 1, 3 and 4, the absence of real roots for every d across five node distributions, and the identity r = Σ L_i p_i relating equations (1) and (4).

One entry of the paper's Table 1 does not reproduce: the sine row at n = 20 is printed as 3.9e−05, but the order 5.5 printed beside it implies 1.7e−2 / 2^5.5 = 3.8e−04, which is what we get, and every other entry in the row matches to two figures. We take it as a misprint.

How it is put together#

• src/matlab/driver.m — reads params.json, calls the method script, writes out.json. It is prepended to whatever is in the editor, so the script's functions become its local functions.
• src/matlab/lib/*.m — the parts that are not the method: node distributions, test functions, the cubic spline, the classical rational interpolant.
• src/methods/*.m — the scripts in the dropdown.
• src/engine/ — the numbl session. One session per script; parameter changes only rewrite params.json and re-run, which is fast enough to drive a slider.
• src/plot/ — a small SVG plotting layer. The palette is documented and was checked with a colour-vision validator rather than by eye.
• src/panels/ — one component per tab.

License#

Apache-2.0, matching numbl.