Inverse quadratic interpolation

Here we look for a root of \(x+\cos(10x)\) that is close to 1.

using FundamentalsNumericalComputation
f = x -> x + cos(10*x);
interval = [0.5,1.5];

plot(f,interval...,grid=:y,label="Function",legend=:bottomright,xlabel="x",ylabel="y")
r = nlsolve(x->f(x[1]),[1.]).zero
1-element Array{Float64,1}:
 0.9678884021293104

We choose three values to get the iteration started.

x = [0.8,1.2,1]
y = @. f(x)
scatter!(x,y,label="initial points")

If we were using “forward” interpolation, we would ask for the polynomial interpolant of \(y\) as a function of \(x\). But that parabola has no real roots.

using Polynomials
q = fit(x,y,2);      # interpolating polynomial
plot!(x->q(x),interval...,l=:dash,label="interpolant",ylim=[-.1,3])

To do inverse interpolation, we swap the roles of \(x\) and \(y\) in the interpolation.

plot(f,interval...,grid=:y,label="Function",legend=:bottomright,xlabel="x",ylabel="y")
scatter!(x,y,label="initial points")

q = fit(y,x,2);      # interpolating polynomial
plot!(y->q(y),y->y,-.1,2.6,l=:dash,label="inverse interpolant")

We seek the value of \(x\) that makes \(y\) zero. This means evaluating \(q\) at zero.

@show x = [x; q(0)];
@show y = [y; f(x[end])];
x = [x; q(0)] = [0.8, 1.2, 1.0, 1.103981385440472]
y = [y; f(x[end])] = [0.6544999661913865, 2.043853958732492, 0.16092847092354756, 1.1482065231940184]

We repeat the process a few more times.

for k = 4:8
    q = fit(y[k-2:k],x[k-2:k],2)
    x = [x; q(0)]
    y = [y; f(x[k+1])]
end

Here is the sequence of errors.

err = @. x - r
9-element Array{Float64,1}:
 -0.16788840212931033
  0.23211159787068958
  0.03211159787068962
  0.1360929833111617
  0.015347343251992385
  0.0032683144541431064
  0.0004617433331599585
  6.295567151926029e-6
  3.1584921345029215e-9

The error seems to be superlinear, but subquadratic.

logerr = @. log(abs(err));
ratios = @. logerr[2:end] / logerr[1:end-1]
8-element Array{Float64,1}:
 0.8184775454303989
 2.3542970942956387
 0.5800188686031708
 2.094252632650976
 1.3702986042756682
 1.341928370566436
 1.5592295620791041
 1.6344120701878158