Notebooks · posts/runge-phenomenon

The Runge phenomenon, where more points give a worse fit

Adam Zettel · September 21, 2026 · written with an agent

Interpolate $1/ (1 + 25 x^{2})$ by a polynomial through equally spaced points and the error grows as you add points. Move the same number of points to $cos (j π / n)$ and it falls geometrically, at a rate you can work out from where the function's poles are.

f (x) = \frac{1}{1 + 25 x ^{2}}, ∥ f - p_{n} ∥_{\infty} \to \infty

The function has poles at $x = \pm i /5$ . Polynomial interpolation in Chebyshev points converges like $ρ^{- n}$ , where $ρ$ is the size of the largest ellipse with foci $\pm 1$ that the function is analytic inside. The ellipse through $\pm i /5$ has $ρ = \frac{1}{5} + 1 + \frac{1}{25} = 1.2198$ . Both interpolants use the barycentric formula, which is stable in Chebyshev points. With equally spaced points the problem is ill-conditioned (its Lebesgue constant grows exponentially with $n$ ), and this formula adds rounding error of its own.

Look inside this program

Atlas
The program’s map: every function and which calls which.
There’s no recorded run of this program yet.

Shift+Enter in a cell runs it and moves to the next. Running a cell runs every cell above it too, as one program.

[1]lines 20–121

pi = 3.141592653589793
M = 400
XS = linspace(-1.0, 1.0, M)

#| Runge's function, on an array of points.
def runge(x) {
    return 1.0 / (1.0 + 25.0 * x * x)
}

#| The degree-n interpolant through the given nodes, by the barycentric formula, on XS.
def interpolate(nodes, weights) {
    n = len(nodes)
    values = runge(nodes)
    out = zeros(M)
    for i to M {
        top = 0.0
        bottom = 0.0
        hit = -1.0
        for j to n {
            d = XS[i] - nodes[j]
            if abs(d) < 1.0e-14 {
                hit = j
            }
            else {
                top = top + weights[j] * values[j] / d
                bottom = bottom + weights[j] / d
            }
        }
        if hit >= 0.0 {
            out[i] = values[hit]
        }
        else {
            out[i] = top / bottom
        }
    }
    return out
}

#| n + 1 equally spaced points, with barycentric weights (-1)^j C(n, j).
def equispaced(n) {
    nodes = linspace(-1.0, 1.0, n + 1.0)
    w = zeros(n + 1.0)
    w[0] = 1.0
    for j to n {
        w[j + 1.0] = -w[j] * (n - j) / (j + 1.0)
    }
    return interpolate(nodes, w)
}

#| n + 1 Chebyshev points cos(j pi / n), with weights (-1)^j, halved at the two ends.
def chebyshev(n) {
    nodes = zeros(n + 1.0)
    w = zeros(n + 1.0)
    sign = 1.0
    for j to n + 1.0 {
        nodes[j] = cos(j * pi / n)
        w[j] = sign
        sign = -sign
    }
    w[0] = 0.5 * w[0]
    w[n] = 0.5 * w[n]
    return interpolate(nodes, w)
}

# a degree-20 interpolant goes through 21 points
DEGREE = 20.0
POINTS = DEGREE + 1.0
# the function's poles are at x = +-i POLE
POLE = 1.0 / 5.0

EXACT = runge(XS)
equi_10 = take the largest difference between equispaced(10.0) and EXACT
equi_20 = take the largest difference between equispaced(DEGREE) and EXACT
cheb_20 = take the largest difference between chebyshev(DEGREE) and EXACT
cheb_40 = take the largest difference between chebyshev(40.0) and EXACT
rate = pow(cheb_20 / cheb_40, 1.0 / 20.0)
predicted_rate = POLE + sqrt(1.0 + POLE * POLE)
show equi_10
show equi_20
show cheb_20
show cheb_40
show rate
show predicted_rate

figure("Interpolating $1/(1+25x^2)$ at 21 points", "x", "interpolant", -1.0, 1.0, -0.4, 1.4)
plot(XS, EXACT, "$1/(1 + 25x^2)$")
plot(XS, equispaced(20.0), "equally spaced")
plot(XS, chebyshev(20.0), "Chebyshev")
style("$1/(1 + 25x^2)$", "dash thick")

test "with equally spaced points, doubling them makes the error more than thirty times larger" {
    assert_gt(equi_20, 30.0 * equi_10)
}
test "with Chebyshev points the same 21 are good to two percent" {
    assert_lt(cheb_20, 0.02)
}
test "and the error falls like rho^-n, with rho set by the poles at +-i/5" {
    assert_near(rate, predicted_rate, 0.02)
}

# SPDX-License-Identifier: MIT OR Apache-2.0
# SPDX-FileCopyrightText: 2026 The Viv authors

pi = 3.141592653589793
M = 400
XS = linspace(-1.0, 1.0, M)

#| Runge's function, on an array of points.
def runge(x) {
    return 1.0 / (1.0 + 25.0 * x * x)
}

#| The degree-n interpolant through the given nodes, by the barycentric formula, on XS.
def interpolate(nodes, weights) {
    n = len(nodes)
    values = runge(nodes)
    out = zeros(M)
    for i to M {
        top = 0.0
        bottom = 0.0
        hit = -1.0
        for j to n {
            d = XS[i] - nodes[j]
            if abs(d) < 1.0e-14 {
                hit = j
            }
            else {
                top = top + weights[j] * values[j] / d
                bottom = bottom + weights[j] / d
            }
        }
        if hit >= 0.0 {
            out[i] = values[hit]
        }
        else {
            out[i] = top / bottom
        }
    }
    return out
}

#| n + 1 equally spaced points, with barycentric weights (-1)^j C(n, j).
def equispaced(n) {
    nodes = linspace(-1.0, 1.0, n + 1.0)
    w = zeros(n + 1.0)
    w[0] = 1.0
    for j to n {
        w[j + 1.0] = -w[j] * (n - j) / (j + 1.0)
    }
    return interpolate(nodes, w)
}

#| n + 1 Chebyshev points cos(j pi / n), with weights (-1)^j, halved at the two ends.
def chebyshev(n) {
    nodes = zeros(n + 1.0)
    w = zeros(n + 1.0)
    sign = 1.0
    for j to n + 1.0 {
        nodes[j] = cos(j * pi / n)
        w[j] = sign
        sign = -sign
    }
    w[0] = 0.5 * w[0]
    w[n] = 0.5 * w[n]
    return interpolate(nodes, w)
}

# a degree-20 interpolant goes through 21 points
DEGREE = 20.0
POINTS = DEGREE + 1.0
# the function's poles are at x = +-i POLE
POLE = 1.0 / 5.0

EXACT = runge(XS)
equi_10 = take the largest difference between equispaced(10.0) and EXACT
equi_20 = take the largest difference between equispaced(DEGREE) and EXACT
cheb_20 = take the largest difference between chebyshev(DEGREE) and EXACT
cheb_40 = take the largest difference between chebyshev(40.0) and EXACT
rate = pow(cheb_20 / cheb_40, 1.0 / 20.0)
predicted_rate = POLE + sqrt(1.0 + POLE * POLE)
show equi_10
show equi_20
show cheb_20
show cheb_40
show rate
show predicted_rate

figure("Interpolating $1/(1+25x^2)$ at 21 points", "x", "interpolant", -1.0, 1.0, -0.4, 1.4)
plot(XS, EXACT, "$1/(1 + 25x^2)$")
plot(XS, equispaced(20.0), "equally spaced")
plot(XS, chebyshev(20.0), "Chebyshev")
style("$1/(1 + 25x^2)$", "dash thick")

test "with equally spaced points, doubling them makes the error more than thirty times larger" {
    assert_gt(equi_20, 30.0 * equi_10)
}
test "with Chebyshev points the same 21 are good to two percent" {
    assert_lt(cheb_20, 0.02)
}
test "and the error falls like rho^-n, with rho set by the poles at +-i/5" {
    assert_near(rate, predicted_rate, 0.02)
}

# SPDX-License-Identifier: MIT OR Apache-2.0
# SPDX-FileCopyrightText: 2026 The Viv authors