Chebyshev rational approximation with Remez algorithm: how do I tell how many digits are accurate?

Question

I have implemented the Remez algorithm for Chebyshev rational approximation using a multiple precision numerics library (mpmath via SymPy). Wikipedia has suggested using $\max\{|z_i|\} - \min\{|z_i|\}$ as a stopping criterion. Here $z_i$ are the values of the error function $E(x) = f(x) - r_i(x)$ evaluated at the local extrema, that is, $z_i = E(x_i)$. Other sources have suggested stopping when the alternates $x_i$ stop changing, but Wikipedia's suggestion seems more reasonable to me, as the theory says that the perfect minimax rational approximation has equioscillating error, and $\max\{|z_i|\} - \min\{|z_i|\}$ is a measure of closeness to that.

Empirically, when I perform the computations with $N$ digits of precision, the convergence happens when $\max\{|z_i|\} - \min\{|z_i|\}$ goes under $10^N$ (that is, iterating past that doesn't produce any further improvement), so I've set that as my stopping condition. If anyone knows why this might be the case (or might not actually be the case) that would be great to know.

My question is, once the approximation has converged using $N$ digits of precision and $\max\{|z_i|\} - \min\{|z_i|\} < \epsilon$, how many digits of the coefficients in final approximation $r(x)$ are correct? Obviously my end goal is to get at least enough digits for some double floating point calculations.

My only goal post thus far has been a paper with published coefficients.

My question is twofold:

1. How can I tell how many digits of the coefficients in my computed approximation are accurate?

2. How can I put a lower bound on the number of accurate digits (preferably by increasing the number of digits of precision in the algorithm, $N$, or by changing the convergence conditions)?

I only care about approximating $e^{-x}$ on $[0, \infty)$, although general answers would be nice. I've been using the technique of the above paper (translate to $[-1, 1]$ via $z = c\frac{1 + x}{1 - x}$ with suitably chosen $c$ to spread out the errors, then translate back with $x = \frac{z - c}{z + c}$).

So far, I've thought of the following:

Suppose we have a polynomial approximation instead of a rational approximation (i.e., $E(x) = f(x) - p(x) < \epsilon$ is the error for the best minimax approximation). If we perturb a coefficient of $p(x)$ by $\delta$, and let $$E^* = f(x) - p^*(x),$$ where $$p^*(x) = a_0 + \cdots + (a_k + \delta)x^k + \cdots + a_nx^n,$$ then $$|E^*| = |E - \delta x^k| <= \epsilon + \delta$$ because $x\in [-1, 1].$ But I'm not sure if this helps me. I also don't know how to apply it to the rational case (how does the denominator affect things?), and finally, I don't know how the inverse transformation $\frac{z - c}{z + c}$ will affect things.

Answer 1

For practical purposes you can compute the approximation for some large $K >> N$, say $K = N + 100$ (I've used $K=1000$). Then, although you do not know how many digits of the $K$-digit approximation are accurate, you can be sure that at least $N$ of them are, so you can then compare how many digits of the $N$-digit approximation agree with the first $N$ digits of the $K$ approximation.

This is a little hand-wavy, because it relies on knowing that the algorithm has properly "converged". I've actually had to adjust the $10^N$ condition I said above to match the published digits in that paper, so this does matter. However, (as far as I can tell) even partial convergence should have the first digits correct, so if $K$ is large enough, then this shouldn't be an issue.