Convergence of AAA rational approximations

using RationalFunctionApproximation, CairoMakie

For a function that is analytic on its domain, the AAA algorithm typically converges at least root-exponentially. In order to observe the convergence, we can construct an approximation that preserves the history of its construction. For example, we approximate an oscillatory function over $[-1,1]$ via:

f = x -> cos(11x)
r = approximate(f, unit_interval, stats=true)
convergenceplot(r)
Example block output

In the plot above, the markers show the estimated max-norm error of the $n$-point AAA rational interpolant over the domain as a function of the numerator/denominator degree $n-1$. The red circles indicate that a pole of the rational interpolant lies on the interval itself. We can verify this by using rewind to recover the degree-7 approximation that was found along the way to r:

julia> r7 = rewind(r, 7)Barycentric rational function of type (7,7) on the domain: Path{Float64} with 1 curve
julia> poles(r7)7-element Vector{ComplexF64}: -0.7865385850382732 - 0.26769746910591097im -0.7865385850382732 + 0.2676974691059109im -0.3174146424804872 + 0.28894380984739326im -0.3174146424804871 - 0.2889438098473932im 0.16129605360549806 + 0.2496796849023167im 0.1612960536054981 - 0.24967968490231668im 0.7365067113018727 + 0.0im

When the AAA iteration encounters an approximation with such undesired poles, or having less accuracy than a predecessor, the AAA iteration simply disregards that approximation and continues–-unless there have been more than a designated number of consecutive failures, at which the best interpolant ever encountered is returned. That interpolant is indicated by the gold halo in the convergence plot above.

When a singularity is very close to the approximation domain, it can cause stagnation and a large number of bad-pole failures:

f = x -> tanh(3000*(x - 1//5))
r = approximate(f, unit_interval, stats=true)
convergenceplot(r)
Example block output

This effect is thought to be mainly due to roundoff and conditioning of the problem. If we use more accurate floating-point arithmetic, we can see that the AAA convergence continues steadily past the previous plateau. In the following, we apply Double64 arithmetic, having used exact rational numbers already in the definition of f:

using DoubleFloats, ComplexRegions
r = approximate(f, Segment{Double64}(-1, 1), stats=true)
convergenceplot(r)
Example block output

In the extreme case of a function with a singularity on the domain, the convergence can be substantially affected:

f = x -> abs(x - 1/8)
r = approximate(f, unit_interval, stats=true)
convergenceplot(r)
Example block output

In such a case, we might get improvement by increasing the number of allowed consecutive failures via the lookahead keyword argument:

r = approximate(f, unit_interval, stats=true, lookahead=20)
convergenceplot(r)
Example block output