2607.19283v1 / ENOTV/Euler.lean

all files

import ENOTV.Compact
import Mathlib.NumberTheory.BernoulliPolynomials

/-!
# Euler polynomials and the matching identity

The even-order counterexample uses the polynomial `E_d` characterized by
`E_d(x+1)+E_d(x)=2x^d`.  We construct it explicitly from Bernoulli
polynomials, so no existence postulate is hidden in the construction.
-/

noncomputable section

open scoped BigOperators
open Polynomial

namespace ENOTV

/-- The Euler polynomial over `ℚ`, in its Bernoulli-polynomial formula. -/
def eulerPoly (d : ℕ) : ℚ[X] :=
  C ((2 : ℚ) ^ (d + 1) / (d + 1)) *
    ((Polynomial.bernoulli (d + 1)).comp (C (1 / 2 : ℚ) * (X + 1)) -
     (Polynomial.bernoulli (d + 1)).comp (C (1 / 2 : ℚ) * X))

theorem eulerPoly_eval (d : ℕ) (x : ℚ) :
    (eulerPoly d).eval x =
      (2 : ℚ) ^ (d + 1) / (d + 1) *
        ((Polynomial.bernoulli (d + 1)).eval ((x + 1) / 2) -
         (Polynomial.bernoulli (d + 1)).eval (x / 2)) := by
  simp only [eulerPoly, eval_mul, eval_C, eval_sub, eval_comp, eval_add, eval_X, eval_one]
  ring

/-- Defining Euler identity. -/
theorem eulerPoly_add_one (d : ℕ) (x : ℚ) :
    (eulerPoly d).eval (x + 1) + (eulerPoly d).eval x = 2 * x ^ d := by
  rw [eulerPoly_eval, eulerPoly_eval]
  have hB := Polynomial.bernoulli_eval_one_add (d + 1) (x / 2)
  rw [show (x + 1 + 1) / 2 = 1 + x / 2 by ring, hB]
  have hd1 : (d + 1 : ℚ) ≠ 0 := by positivity
  rw [Nat.add_sub_cancel, div_pow]
  field_simp [hd1]
  push_cast
  ring

theorem eulerPoly_polynomial_identity (d : ℕ) :
    (eulerPoly d).comp (X + 1) + eulerPoly d = C 2 * X ^ d := by
  apply Polynomial.funext
  intro x
  simpa [Polynomial.eval_comp, add_comm] using eulerPoly_add_one d x

theorem eulerPoly_natDegree_le (d : ℕ) : (eulerPoly d).natDegree ≤ d := by
  by_contra hle
  have hdn : d < (eulerPoly d).natDegree := Nat.lt_of_not_ge hle
  let n := (eulerPoly d).natDegree
  have hn0 : n ≠ 0 := by omega
  have hqdeg : (X + (1 : ℚ[X])).natDegree = 1 := by
    simpa only [Polynomial.C_1] using Polynomial.natDegree_X_add_C (1 : ℚ)
  have htrans :
      ((eulerPoly d).comp (X + 1)).natDegree = n := by
    rw [Polynomial.natDegree_comp, hqdeg, mul_one]
  have hlc :
      ((eulerPoly d).comp (X + 1)).leadingCoeff =
        (eulerPoly d).leadingCoeff := by
    rw [Polynomial.leadingCoeff_comp (by omega), show X + (1 : ℚ[X]) = X + C 1 by simp,
      Polynomial.leadingCoeff_X_add_C, one_pow, mul_one]
  have hcoeffcomp :
      ((eulerPoly d).comp (X + 1)).coeff n = (eulerPoly d).leadingCoeff := by
    rw [← htrans, Polynomial.coeff_natDegree, hlc]
  have hcoeffself :
      (eulerPoly d).coeff n = (eulerPoly d).leadingCoeff := by
    exact Polynomial.coeff_natDegree
  have hc := congrArg (fun p : ℚ[X] => p.coeff n)
    (eulerPoly_polynomial_identity d)
  simp only [Polynomial.coeff_add] at hc
  rw [hcoeffcomp, hcoeffself] at hc
  have hnD : n ≠ d := by omega
  rw [Polynomial.coeff_C_mul_X_pow] at hc
  simp [hnD] at hc
  have hp : eulerPoly d ≠ 0 := by
    intro hp
    simp [n, hp] at hn0
  exact hp hc

/-- For even degree the Euler polynomial is symmetric about `1/2`. -/
theorem eulerPoly_one_sub (d : ℕ) (hd : Even d) (x : ℚ) :
    (eulerPoly d).eval (1 - x) = (eulerPoly d).eval x := by
  have hodd : Odd (d + 1) := hd.add_one
  have h₁ := Polynomial.bernoulli_eval_one_sub (d + 1) (x / 2)
  have h₂ := Polynomial.bernoulli_eval_one_sub (d + 1) ((x + 1) / 2)
  rw [hodd.neg_one_pow] at h₁ h₂
  rw [eulerPoly_eval, eulerPoly_eval]
  rw [show (1 - x + 1) / 2 = 1 - x / 2 by ring, h₁,
    show (1 - x) / 2 = 1 - (x + 1) / 2 by ring, h₂]
  ring

theorem eulerPoly_zero (d : ℕ) (hd : 0 < d) (heven : Even d) :
    (eulerPoly d).eval 0 = 0 := by
  have hsym := eulerPoly_one_sub d heven 0
  have hid := eulerPoly_add_one d 0
  norm_num [zero_pow hd.ne'] at hid
  norm_num at hsym
  linarith

theorem eulerPoly_one (d : ℕ) (hd : 0 < d) (heven : Even d) :
    (eulerPoly d).eval 1 = 0 := by
  simpa using (eulerPoly_one_sub d heven 0).trans (eulerPoly_zero d hd heven)

/-- The rescaled block polynomial `P_L(t)=(2L)^d E_d(t/(2L))`. -/
def blockPoly (d L : ℕ) (t : ℚ) : ℚ :=
  (2 * L : ℚ) ^ d * (eulerPoly d).eval (t / (2 * L))

/-- Polynomial whose evaluation is `blockPoly`. -/
def blockPolynomial (d L : ℕ) : ℚ[X] :=
  C ((2 * L : ℚ) ^ d) *
    (eulerPoly d).comp (C ((2 * L : ℚ)⁻¹) * X)

theorem blockPolynomial_eval (d L : ℕ) (t : ℚ) :
    (blockPolynomial d L).eval t = blockPoly d L t := by
  by_cases hL : L = 0
  · subst L
    simp [blockPolynomial, blockPoly]
  · simp [blockPolynomial, blockPoly]
    congr 2
    field_simp <;> simp [hL]

theorem blockPolynomial_natDegree_le (d L : ℕ) :
    (blockPolynomial d L).natDegree ≤ d := by
  calc
    (blockPolynomial d L).natDegree ≤
        ((eulerPoly d).comp (C ((2 * L : ℚ)⁻¹) * X)).natDegree := by
      exact Polynomial.natDegree_C_mul_le _ _
    _ = (eulerPoly d).natDegree *
        (C ((2 * L : ℚ)⁻¹) * X).natDegree := Polynomial.natDegree_comp
    _ ≤ (eulerPoly d).natDegree * 1 := by
      gcongr
      simpa only [pow_one] using
        Polynomial.natDegree_C_mul_X_pow_le ((2 * L : ℚ)⁻¹) 1
    _ ≤ d := by simpa using eulerPoly_natDegree_le d

theorem blockPolynomial_high_difference (d L : ℕ) :
    (fwdDiff (1 : ℚ))^[d + 1] (blockPolynomial d L).eval = 0 := by
  apply Polynomial.fwdDiff_iter_eq_zero_of_degree_lt
  exact (blockPolynomial_natDegree_le d L).trans_lt (Nat.lt_succ_self d)

theorem blockPoly_zero (d L : ℕ) (hd : 0 < d) (heven : Even d) :
    blockPoly d L 0 = 0 := by
  simp [blockPoly, eulerPoly_zero d hd heven]

theorem blockPoly_two_mul (d L : ℕ) (hd : 0 < d) (heven : Even d) :
    blockPoly d L (2 * L) = 0 := by
  by_cases hL : L = 0
  · simp [hL, blockPoly, zero_pow hd.ne']
  · simp [blockPoly, hL, eulerPoly_one d hd heven]

/-- Exact, scale-independent matching at an alternating-block junction. -/
theorem blockPoly_matching (d L : ℕ) (hL : 0 < L) (h : ℚ) :
    blockPoly d L (2 * L + h) + blockPoly d L h = 2 * h ^ d := by
  rw [blockPoly, blockPoly]
  have hid := eulerPoly_add_one d (h / (2 * L))
  have hden : (2 * (L : ℚ)) ≠ 0 := by positivity
  rw [show (2 * (L : ℚ) + h) / (2 * (L : ℚ)) = h / (2 * L) + 1 by
      field_simp
      <;> ring]
  calc
    (2 * (L : ℚ)) ^ d * (eulerPoly d).eval (h / (2 * L) + 1) +
          (2 * (L : ℚ)) ^ d * (eulerPoly d).eval (h / (2 * L)) =
        (2 * (L : ℚ)) ^ d *
          ((eulerPoly d).eval (h / (2 * L) + 1) +
           (eulerPoly d).eval (h / (2 * L))) := by ring
    _ = (2 * (L : ℚ)) ^ d * (2 * (h / (2 * L)) ^ d) := by rw [hid]
    _ = 2 * h ^ d := by
      rw [div_pow]
      field_simp [hden]

end ENOTV