Optimal polynomial approximation of rational function $\frac{1}{1-x}$

Question

I've been working on the following polynomial approximation problem. I want to find the optimal Chebyshev approximation of the rational function $\frac{1}{1-x}$ on the real interval $x\in[-\rho, \rho]$, $\rho<1$ using $T$-degree polynomials $h_T(x)$. I formulate the approximation problem as the following optimization problem: $$ \min_{\alpha_i} \left\|\frac{1}{1-x}-\alpha_Tx^T-\alpha_{T-1}x^{T-1}-\cdots-\alpha_1x-\alpha_0\right\|_{\infty, \rho}, $$ where the notation $\lVert\cdot\rVert_{\infty, \rho}$ stands for the Chebyshev norm of polynomials defined as: $$ \lVert f(x)\rVert_{\infty, \rho}=\sup_{x\in[-\rho, \rho]} \lvert f(x)\rvert. $$ Currently, I can solve the above problem numerically using the Pólya–Szegő theorem. To be specific, the above optimization problem can be massaged as: $$ \begin{aligned} & \min_{\alpha_i}\ \delta, \\ \text{s.t. } & \frac{1}{1-x}-\sum_{i=0}^T\alpha_ix^i\leq \delta, \forall x\in[-\rho, \rho], \\ & \frac{1}{1-x}-\sum_{i=0}^T \alpha_ix^i\geq -\delta, \forall x\in[-\rho, \rho], \end{aligned} $$ which could be further modified to: $$ \begin{aligned} & \min_{\alpha_i}\ \delta, \\ s.t.\ & -1+(1-x)\left(\sum_{i=0}^T\alpha_ix^i+\delta\right)\geq 0, \forall x\in[-\rho, \rho], \\ & 1+(1-x)\left(-\sum_{i=0}^T \alpha_ix^i+\delta\right)\geq 0, \forall x\in[-\rho, \rho]. \end{aligned} $$ According to Powers and Reznick - Polynomials that are positive on an interval, the inequalities of the optimization problem hold as long as the polynomials on the LHS can be written as $(f(x))^2+(\rho^2-x^2)(g(x))^2$, where $f(x)$ is a $T$-degree polynomial and $g(x)$ is a ($T-1$)-degree polynomial.

Using this technique, I massaged the above optimization problem in the polynomial form to an equivalent Semi-Definite Programming (SDP) problem and solved the problem using SDP solvers.

Here are some numerical results I've got so far. First, I observe that the approximation error decreases exponentially w.r.t. the degree of the polynomial $h_T(x)$, and my simulation result is It can be seen that the approximation error attenuates exponentially w.r.t. $T$ before $T=12$. I believe that my SDP solver may suffer from numerical issues when the approximation error $\epsilon\leq 1e-7$ (maybe because I set the duality gap of the solver to be $1e-6$) or when the dimension of the matrices are large so that the curve looks weird after $T=12$.

Moreover, I study the relationship between the approximation error and the length of the interval $\rho$ plotted here (the degree of $h_T(x)$ is fixed to be $T=6$). It seems that the relationship is approximately (but not exactly) exponential.

On the other hand, I am also interested in the magnitude of the coefficients $\lVert\boldsymbol{\alpha}\rVert_1$, where $\boldsymbol{\alpha}=\left[\alpha_0, \alpha_1, \dotsc, \alpha_T\right]$. The magnitude $\lVert\boldsymbol{\alpha}\rVert_1$ w.r.t. $T$ are visualized in and $\lVert\boldsymbol{\alpha}\rVert_1$ w.r.t. $\rho$ here (still $T=6$). It can be seen that $\lVert\boldsymbol{\alpha}\rVert_1$ grows linearly with $T$ for a small $\rho$, and grows exponentially w.r.t. $T$ for a large $\rho$ ($\rho\geq 0.7$). Besides, $\lVert\boldsymbol{\alpha}\rVert_1$ increases quicker than exponentially w.r.t. $\rho$.

Despite the numerical results so far, the problem lacks theoretical guarantee on both the approximation error $\epsilon$ and the magnitude of the coefficients $\lVert\boldsymbol{\alpha}\rVert_1$. Could anyone help me to find a way out?

Answer 1

This problem has an exact solution, written in the book

N. I. Akhiezer, Theory of approximation. Dover Publications, Inc., New York, 1992, Chap II section 37. The error is $$\frac{(1-\sqrt{1-\rho^2})^n}{\rho^{n-2}(1-\rho^2)}\sim 2^{-n}\rho^{n+2}/(1-\rho^2),\quad n\to\infty.$$ The LHS is the exact expression of the error. The polynomial of best approximation is unique and is also described, though the formula is too complicated to be reproduced here.

Remark. Russian original of Akhiezer's book is available on my web page: https://www.math.purdue.edu/~eremenko/books-papers.html (An English translation exists.). The result is on p. 70-71, and the polynomial of the best approximation is written there.

Answer 2

Your problem is one of approximation theory (as one of your tags states). That wikipedia page has some links to some good general techniques, namely the Remez algorithm for provably optimal $\epsilon$ (though no guarantees regarding $\lVert \alpha\rVert_1$), and expansion in Chebyshev series, where you should be quite competitive for both $\epsilon$, and you should be able to get theoretical guarantees for both $\epsilon$ and $\lVert \alpha\rVert_1$.

Specific to (something related to) your function is the concept of division algorithms. These assume integer inputs $N, D$, and approximates the quotient $N/D$. Two common techniques are

Newton-Raphson:

• After some preprocessing to choose an initial estimate $X_0$, compute $X_{i+1} = X_i(2-DX_i)$. Typically $X_3$ or $X_4$ suffices for reasonable $\epsilon$ --- these should be degree $2^3$ or $2^4$ polynomials in $X_0$.

Goldschmidt:

Roughly speaking, for $x\in(1/2, 1]$ one can compute $\frac{N}{D} = \frac{N (1+x)(1+x^2)\dots(1+x^{2^{(n-1)}})}{1-x^{2^n}}\approx N\prod_{i=0}^{n-1}(1+x^{2^i}).$ Again it should suffice for $n$ to be relatively small (though I don't know your precise requirements).

For both of these techniques provable bounds on $\epsilon$ are known, and it should be straightforward to compute $\lVert \alpha\rVert_1$ (or whatever you need).

It's worth mentioning sometimes when additional constraints are put in place, there are more efficient things you can do. Your constraint on $\lVert \alpha\rVert_1$ is vaguely similar to what appears when one wants to give polynomial approximations for homomorphic evaluation using Fully Homomorphic Encryption. See for example this paper, though for the function $\mathsf{sgn}(x)$. That being said, I don't think people in this area have any work that is division-focused (and instead I think Goldschmidt's is considered the best thing to do, but perhaps the bigger takeaway is "don't divide in this computational paradigm", so saying anything is "the best" is perhaps strong).