Why do these two methods give different answers for this simple approximation?

Question

I tried to find Padé approximant of the function below using two different methods, but the results were not equal.

$f(x)=1+x+x^2+x^3+\cdots$

Method 1: Using the direct built in function of Padé approximant.

F = 1 + x + x^2 + x^3 + x^4 + x^5 + x^6;
PadeApproximant[F, {x, 0, {3, 2}}]

The output will be $1/(1 - x)$

Method 2: Using Eq .15 on https://mathworld.wolfram.com/PadeApproximant.html. Here P is the numerator of Eq .15, and it becomes zero.

P = Det[{{1, 1, 1}, {1, 1, 1}, {x^2 + x^3, x^1 + x^2 + x^3, 1 + x + x^2 + x^3}}]

I have tried it for some functions such as $e^x$ and I got correct equal answers from both methods. I do not know to fix the second method for $f(x)$. Thank you in advance.

Answer 1

To summarize, the two results are different because the second method uses equations that are linearly dependent for some functions. By themselves, they do not contain enough information to determine a unique solution.

Long answer: Take a look at Eqs 12-14 in the link that you cited and imagine all of the $a_i=1$, as in your example. All of those equations would then be the same, not independent. Now look at Eq 15 and imagine all of the $a_i=1$. The first $M$ rows would be identical, so the determinant would be zero.

But, it's even worse than that. You don't have to imagine all of the $a_i=1$. It is only $a_{L-M+1} \dots a_{L+M}$ that appear in Eqs 12-14. So in your example, you could change $a_0$ and $a_1$ and get the same non-answer. And, of course, it's the $a_i$'s being equal to one another, not necessarily equal to 1, that makes the determinant zero.

The following code generates Eqs 7-11 and Eqs 12-14 for the subject polynomial function and attempts to solve the system of equations.

l = 3; m = 2;
rules = {q[0] -> 1, a[_?(# < 0 &)] -> 0, q[_?(# > m &)] -> 0};
eqns7$11 = Table[Sum[a[j] q[irow - j], {j, 0, irow}] == p[irow],
    {irow, 0, l}] /. rules;
eqns12$14 = Table[Sum[a[l - j + irow] q[j], {j, 0, m}] == 0,
    {irow, m}] /. rules;

(eqns = Join[eqns7$11, eqns12$14] /. a[_] -> 1) // Column
soln = First@Solve[eqns]

The Solve solution allows either $q_1$ or $q_2$ to be a free parameter. Setting $q_2=0$ gives $q_1=-1$, in agreement with the PadeApproximant result.

Let's try different choices of q[1] and see what we get for the approximant

    ϕ = Sum[p[k] x^k, {k, 0, l}]/Sum[q[k] x^k, {k, 0, m}]
    ϕ //. Flatten@{q[0] -> 1, soln, q[1] -> -1}
    ϕ //. Flatten@{q[0] -> 1, soln, q[1] -> 0} // Simplify
    ϕ //. Flatten@{q[0] -> 1, soln, q[1] -> 1} // Simplify

(* 
(p[0] + x p[1] + x^2 p[2] + x^3 p[3])/(q[0] + x q[1] + x^2 q[2])
1/(1 - x)
1/(1 - x)
1/(1 - x)  *)

Simplify was used in the above to cancel common polynomial factors in numerator and the denominator.