hex

8.6. Padé approximation🔗

Let s be a power series known to N = m + n + 1 terms. An [m/n] Padé approximant of s is a pair of polynomials p of degree at most m and q ≠ 0 of degree at most n with q s ≡ p modulo xᴺ. One always exists: the congruence is N linear equations in the N + 1 coefficients of p and q. The extended Euclidean algorithm applied to xᴺ and s finds one. Each remainder it produces has the form u xᴺ + v s, so v s is congruent to the remainder modulo xᴺ, and the first remainder of degree at most m gives p and q = v. Hex.DensePoly.padeHomogeneous stops the half-gcd algorithm at that remainder, and returns a Hex.DensePoly.PadeApproximant, which carries proofs of the degree bounds, of q ≠ 0 and of the congruence:

🔗structure
Hex.DensePoly.PadeApproximant.{u} {F : Type u} [DecidableEq F] [Lean.Grind.Field F] {k : ℕ} (s : Hex.TSeries F k) (m n : ℕ) : Type u
Hex.DensePoly.PadeApproximant.{u} {F : Type u} [DecidableEq F] [Lean.Grind.Field F] {k : ℕ} (s : Hex.TSeries F k) (m n : ℕ) : Type u

A homogeneous Padé approximant of numerator degree at most m and denominator degree at most n. A nonzero denominator is stronger than the nontrivial-pair condition and is what the half-gcd construction supplies.

p : Hex.DensePoly F

Numerator polynomial.

q : Hex.DensePoly F

Denominator polynomial.

p_size : self.p.size ≤ m + 1

The numerator has degree at most the requested bound.

q_size : self.q.size ≤ n + 1

The denominator has degree at most the requested bound.

q_ne : self.q ≠ 0

The denominator is nonzero.

congruent : ∀ i < m + n + 1, (self.q * Hex.DensePoly.polyOfSeries s - self.p).coeff i = 0

The denominator times the input agrees with the numerator to the required order.

The pair need not have q(0) ≠ 0. When it does, q is invertible as a power series and s ≡ p / q modulo xᴺ; this is the approximation of the introduction. Hex.DensePoly.pade? scales the denominator to have constant term 1, as in the introduction, and returns a Hex.DensePoly.NormalizedPade:

🔗theorem
Hex.DensePoly.pade?_eq_none_iff.{u} {F : Type u} [DecidableEq F] [Lean.Grind.Field F] (plan : Hex.DensePoly.MulPlan F) {k : ℕ} (s : Hex.TSeries F k) (m n : ℕ) : Hex.DensePoly.pade? plan s m n = none ↔ ¬Nonempty (Hex.DensePoly.NormalizedPade s m n)
Hex.DensePoly.pade?_eq_none_iff.{u} {F : Type u} [DecidableEq F] [Lean.Grind.Field F] (plan : Hex.DensePoly.MulPlan F) {k : ℕ} (s : Hex.TSeries F k) (m n : ℕ) : Hex.DensePoly.pade? plan s m n = none ↔ ¬Nonempty (Hex.DensePoly.NormalizedPade s m n)

Normalized Padé construction fails exactly when no denominator with constant coefficient one satisfies the requested bounds and congruence.

The [2/2] approximant of eˣ is (1 + x/2 + x²/12) / (1 - x/2 + x²/12):

def expSeries : TSeries Rat 5 := TSeries.ofFn fun i => [1, 1, 1/2, 1/6, 1/24].getD i 0 some (#p[1, (1 : Rat)/2, (1 : Rat)/12], #p[1, (-1 : Rat)/2, (1 : Rat)/12])#eval (pade? (karatsubaPlan 32) expSeries 2 2).map fun a => (a.p, a.q)
some (#p[1, (1 : Rat)/2, (1 : Rat)/12], #p[1, (-1 : Rat)/2, (1 : Rat)/12])

For the series x and m = 0, n = 1, the congruence q x ≡ p modulo x² forces p = 0 and q to be a multiple of x. There is no approximant with q(0) = 1, and pade? returns none. The homogeneous approximant is p = 0, q = -x:

def xSeries : TSeries Rat 2 := TSeries.X #guard (pade? (karatsubaPlan 32) xSeries 0 1).isNone #guard (padeHomogeneous (karatsubaPlan 32) xSeries 0 1).q = #p[0, -1] end HexPolyFastChapter