Derive derivatives from Newton forward

Newton's forward interpolation formula, with p = (x − x₀)/h so that dp/dx = 1/h, is

y = y₀ + pΔy₀ + p(p−1)/2!·Δ²y₀ + p(p−1)(p−2)/3!·Δ³y₀ + p(p−1)(p−2)(p−3)/4!·Δ⁴y₀ + …

Expanding the coefficients in powers of p:

y = y₀ + pΔy₀ + (p²−p)/2·Δ²y₀ + (p³−3p²+2p)/6·Δ³y₀ + (p⁴−6p³+11p²−6p)/24·Δ⁴y₀ + …

First derivative

Differentiate with respect to p and use dy/dx = (dy/dp)(dp/dx) = (1/h)(dy/dp):

dy/dx = (1/h)[ Δy₀ + (2p−1)/2·Δ²y₀ + (3p²−6p+2)/6·Δ³y₀ + (4p³−18p²+22p−6)/24·Δ⁴y₀ + … ]

Second derivative

Differentiate again: d²y/dx² = (1/h²)·d²y/dp²:

d²y/dx² = (1/h²)[ Δ²y₀ + (p−1)Δ³y₀ + (6p²−18p+11)/12·Δ⁴y₀ + … ]

At the tabulated point x = x₀ (p = 0)

y′(x₀)  = (1/h)[ Δy₀ − Δ²y₀/2 + Δ³y₀/3 − Δ⁴y₀/4 + … ]
y″(x₀)  = (1/h²)[ Δ²y₀ − Δ³y₀ + 11Δ⁴y₀/12 − … ]

These formulas are used when x is near the beginning of the table; the corresponding backward formulas are used near the end.