Skip to content

Commit

Permalink
feat: define Newton's method and prove decomposition as nilpotent + r…
Browse files Browse the repository at this point in the history
…oot (#10284)

This is just a modified version of the code provided by Antoine Chambert-Loir here: https://leanprover.zulipchat.com/#narrow/stream/116395-maths/topic/jordan-chevalley.20decomposition/near/411402670

Co-authored-by: Antoine Chambert-Loir <antoine.chambert-loir@math.univ-paris-diderot.fr>
  • Loading branch information
ocfnash and Antoine Chambert-Loir committed Feb 6, 2024
1 parent cafe049 commit 229b0ba
Show file tree
Hide file tree
Showing 6 changed files with 199 additions and 21 deletions.
1 change: 1 addition & 0 deletions Mathlib.lean
Expand Up @@ -2096,6 +2096,7 @@ import Mathlib.Dynamics.FixedPoints.Basic
import Mathlib.Dynamics.FixedPoints.Topology
import Mathlib.Dynamics.Flow
import Mathlib.Dynamics.Minimal
import Mathlib.Dynamics.Newton
import Mathlib.Dynamics.OmegaLimit
import Mathlib.Dynamics.PeriodicPts
import Mathlib.FieldTheory.AbelRuffini
Expand Down
128 changes: 128 additions & 0 deletions Mathlib/Dynamics/Newton.lean
@@ -0,0 +1,128 @@
/-
Copyright (c) 2024 Oliver Nash. All rights reserved.
Released under Apache 2.0 license as described in the file LICENSE.
Authors: Antoine Chambert-Loir, Oliver Nash
-/
import Mathlib.Data.Polynomial.AlgebraMap
import Mathlib.Data.Polynomial.Identities
import Mathlib.RingTheory.Nilpotent
import Mathlib.RingTheory.Polynomial.Nilpotent
import Mathlib.RingTheory.Polynomial.Tower

/-!
# Newton-Raphson method
Given a single-variable polynomial `P` with derivative `P'`, Newton's method concerns iteration of
the rational map: `x ↦ x - P(x) / P'(x)`.
Over a field it can serve as a root-finding algorithm. It is also useful tool in certain proofs
such as Hensel's lemma and Jordan-Chevalley decomposition.
## Main definitions / results:
* `Polynomial.newtonMap`: the map `x ↦ x - P(x) / P'(x)`, where `P'` is the derivative of the
polynomial `P`.
* `Polynomial.isFixedPt_newtonMap_of_isUnit_iff`: `x` is a fixed point for Newton iteration iff
it is a root of `P` (provided `P'(x)` is a unit).
* `Polynomial.exists_unique_nilpotent_sub_and_aeval_eq_zero`: if `x` is almost a root of `P` in the
sense that that `P(x)` is nilpotent (and `P'(x)` is a unit) then we may write `x` as a sum
`x = n + r` where `n` is nilpotent and `r` is a root of `P`. This can be used to prove the
Jordan-Chevalley decomposition of linear endomorphims.
-/

open Set Function

noncomputable section

namespace Polynomial

variable {R S : Type*} [CommRing R] [CommRing S] [Algebra R S] (P : R[X]) {x : S}

/-- Given a single-variable polynomial `P` with derivative `P'`, this is the map:
`x ↦ x - P(x) / P'(x)`. When `P'(x)` is not a unit we use a junk-value pattern and send `x ↦ x`. -/
def newtonMap (x : S) : S :=
x - (Ring.inverse <| aeval x (derivative P)) * aeval x P

theorem newtonMap_apply :
P.newtonMap x = x - (Ring.inverse <| aeval x (derivative P)) * (aeval x P) :=
rfl

variable {P}

theorem newtonMap_apply_of_isUnit (h : IsUnit <| aeval x (derivative P)) :
P.newtonMap x = x - h.unit⁻¹ * aeval x P := by
simp [newtonMap_apply, Ring.inverse, h]

theorem newtonMap_apply_of_not_isUnit (h : ¬ (IsUnit <| aeval x (derivative P))) :
P.newtonMap x = x := by
simp [newtonMap_apply, Ring.inverse, h]

theorem isNilpotent_iterate_newtonMap_sub_of_isNilpotent (h : IsNilpotent <| aeval x P) (n : ℕ) :
IsNilpotent <| P.newtonMap^[n] x - x := by
induction n with
| zero => simp
| succ n ih =>
rw [iterate_succ', comp_apply, newtonMap_apply, sub_right_comm]
refine (Commute.all _ _).isNilpotent_sub ih <| (Commute.all _ _).isNilpotent_mul_right ?_
simpa using Commute.isNilpotent_add (Commute.all _ _)
(isNilpotent_aeval_sub_of_isNilpotent_sub P ih) h

theorem isFixedPt_newtonMap_of_aeval_eq_zero (h : aeval x P = 0) :
IsFixedPt P.newtonMap x := by
rw [IsFixedPt, newtonMap_apply, h, mul_zero, sub_zero]

theorem isFixedPt_newtonMap_of_isUnit_iff (h : IsUnit <| aeval x (derivative P)) :
IsFixedPt P.newtonMap x ↔ aeval x P = 0 := by
rw [IsFixedPt, newtonMap_apply, sub_eq_self, Ring.inverse_mul_eq_iff_eq_mul _ _ _ h, mul_zero]

/-- This is really an auxiliary result, en route to
`Polynomial.exists_unique_nilpotent_sub_and_aeval_eq_zero`. -/
theorem aeval_pow_two_pow_dvd_aeval_iterate_newtonMap
(h : IsNilpotent (aeval x P)) (h' : IsUnit (aeval x <| derivative P)) (n : ℕ) :
(aeval x P) ^ (2 ^ n) ∣ aeval (P.newtonMap^[n] x) P := by
induction n with
| zero => simp
| succ n ih =>
have ⟨d, hd⟩ := binomExpansion (P.map (algebraMap R S)) (P.newtonMap^[n] x)
(-Ring.inverse (aeval (P.newtonMap^[n] x) <| derivative P) * aeval (P.newtonMap^[n] x) P)
rw [eval_map_algebraMap, eval_map_algebraMap] at hd
rw [iterate_succ', comp_apply, newtonMap_apply, sub_eq_add_neg, neg_mul_eq_neg_mul, hd]
refine dvd_add ?_ (dvd_mul_of_dvd_right ?_ _)
· convert dvd_zero _
have : IsUnit (aeval (P.newtonMap^[n] x) <| derivative P) :=
isUnit_aeval_of_isUnit_aeval_of_isNilpotent_sub h' <|
isNilpotent_iterate_newtonMap_sub_of_isNilpotent h n
rw [derivative_map, eval_map_algebraMap, ← mul_assoc, mul_neg, Ring.mul_inverse_cancel _ this,
neg_mul, one_mul, add_right_neg]
· rw [neg_mul, even_two.neg_pow, mul_pow, pow_succ', pow_mul]
exact dvd_mul_of_dvd_right (pow_dvd_pow_of_dvd ih 2) _

/-- If `x` is almost a root of `P` in the sense that that `P(x)` is nilpotent (and `P'(x)` is a
unit) then we may write `x` as a sum `x = n + r` where `n` is nilpotent and `r` is a root of `P`.
Moreover, `n` and `r` are unique.
This can be used to prove the Jordan-Chevalley decomposition of linear endomorphims. -/
theorem exists_unique_nilpotent_sub_and_aeval_eq_zero
(h : IsNilpotent (aeval x P)) (h' : IsUnit (aeval x <| derivative P)) :
∃! r, IsNilpotent (x - r) ∧ aeval r P = 0 := by
simp_rw [(neg_sub _ x).symm, isNilpotent_neg_iff]
refine exists_unique_of_exists_of_unique ?_ fun r₁ r₂ ⟨hr₁, hr₁'⟩ ⟨hr₂, hr₂'⟩ ↦ ?_
· -- Existence
obtain ⟨n, hn⟩ := id h
refine ⟨P.newtonMap^[n] x, isNilpotent_iterate_newtonMap_sub_of_isNilpotent h n, ?_⟩
rw [← zero_dvd_iff, ← pow_eq_zero_of_le n.lt_two_pow.le hn]
exact aeval_pow_two_pow_dvd_aeval_iterate_newtonMap h h' n
· -- Uniqueness
have ⟨u, hu⟩ := binomExpansion (P.map (algebraMap R S)) r₁ (r₂ - r₁)
suffices IsUnit (aeval r₁ (derivative P) + u * (r₂ - r₁)) by
rwa [derivative_map, eval_map_algebraMap, eval_map_algebraMap, eval_map_algebraMap,
add_sub_cancel'_right, hr₂', hr₁', zero_add, pow_two, ← mul_assoc, ← add_mul, eq_comm,
this.mul_right_eq_zero, sub_eq_zero, eq_comm] at hu
have : IsUnit (aeval r₁ (derivative P)) :=
isUnit_aeval_of_isUnit_aeval_of_isNilpotent_sub h' hr₁
rw [← sub_sub_sub_cancel_right r₂ r₁ x]
refine IsNilpotent.isUnit_add_left_of_commute ?_ this (Commute.all _ _)
exact (Commute.all _ _).isNilpotent_mul_right <| (Commute.all _ _).isNilpotent_sub hr₂ hr₁

end Polynomial
53 changes: 38 additions & 15 deletions Mathlib/RingTheory/Nilpotent.lean
Expand Up @@ -109,22 +109,45 @@ lemma IsNilpotent.map_iff [MonoidWithZero R] [MonoidWithZero S] {r : R} {F : Typ
IsNilpotent (f r) ↔ IsNilpotent r :=
fun ⟨k, hk⟩ ↦ ⟨k, (map_eq_zero_iff f hf).mp <| by rwa [map_pow]⟩, fun h ↦ h.map f⟩

theorem IsNilpotent.sub_one_isUnit [Ring R] {r : R} (hnil : IsNilpotent r) : IsUnit (r - 1) := by
obtain ⟨n, hn⟩ := hnil
refine' ⟨⟨r - 1, -∑ i in Finset.range n, r ^ i, _, _⟩, rfl⟩
· rw [mul_neg, mul_geom_sum, hn]
simp
· rw [neg_mul, geom_sum_mul, hn]
simp

theorem Commute.IsNilpotent.add_isUnit [Ring R] {r : R} {u : Rˣ} (hnil : IsNilpotent r)
(hru : Commute r (↑u⁻¹ : R)) : IsUnit (u + r) := by
rw [← Units.isUnit_mul_units _ u⁻¹, add_mul, Units.mul_inv, ← IsUnit.neg_iff, add_comm, neg_add,
← sub_eq_add_neg]
theorem IsUnit.isNilpotent_mul_unit_of_commute_iff [MonoidWithZero R] {r u : R}
(hu : IsUnit u) (h_comm : Commute r u) :
IsNilpotent (r * u) ↔ IsNilpotent r :=
exists_congr fun n ↦ by rw [h_comm.mul_pow, (hu.pow n).mul_left_eq_zero]

theorem IsUnit.isNilpotent_unit_mul_of_commute_iff [MonoidWithZero R] {r u : R}
(hu : IsUnit u) (h_comm : Commute r u) :
IsNilpotent (u * r) ↔ IsNilpotent r :=
h_comm ▸ hu.isNilpotent_mul_unit_of_commute_iff h_comm

theorem IsNilpotent.isUnit_sub_one [Ring R] {r : R} (hnil : IsNilpotent r) : IsUnit (r - 1) := by
obtain ⟨n, hn⟩ := hnil
refine' IsNilpotent.sub_one_isUnit ⟨n, _⟩
rw [neg_pow, hru.mul_pow, hn]
simp
refine ⟨⟨r - 1, -∑ i in Finset.range n, r ^ i, ?_, ?_⟩, rfl⟩
· simp [mul_geom_sum, hn]
· simp [geom_sum_mul, hn]

theorem IsNilpotent.isUnit_one_sub [Ring R] {r : R} (hnil : IsNilpotent r) : IsUnit (1 - r) := by
rw [← IsUnit.neg_iff, neg_sub]
exact isUnit_sub_one hnil

theorem IsNilpotent.isUnit_add_one [Ring R] {r : R} (hnil : IsNilpotent r) : IsUnit (r + 1) := by
rw [← IsUnit.neg_iff, neg_add']
exact isUnit_sub_one hnil.neg

theorem IsNilpotent.isUnit_one_add [Ring R] {r : R} (hnil : IsNilpotent r) : IsUnit (1 + r) :=
add_comm r 1 ▸ isUnit_add_one hnil

theorem IsNilpotent.isUnit_add_left_of_commute [Ring R] {r u : R}
(hnil : IsNilpotent r) (hu : IsUnit u) (h_comm : Commute r u) :
IsUnit (u + r) := by
rw [← Units.isUnit_mul_units _ hu.unit⁻¹, add_mul, IsUnit.mul_val_inv]
replace h_comm : Commute r (↑hu.unit⁻¹) := Commute.units_inv_right h_comm
refine IsNilpotent.isUnit_one_add ?_
exact (hu.unit⁻¹.isUnit.isNilpotent_mul_unit_of_commute_iff h_comm).mpr hnil

theorem IsNilpotent.isUnit_add_right_of_commute [Ring R] {r u : R}
(hnil : IsNilpotent r) (hu : IsUnit u) (h_comm : Commute r u) :
IsUnit (r + u) :=
add_comm r u ▸ hnil.isUnit_add_left_of_commute hu h_comm

section NilpotencyClass

Expand Down
2 changes: 1 addition & 1 deletion Mathlib/RingTheory/Polynomial/IrreducibleRing.lean
Expand Up @@ -57,5 +57,5 @@ theorem Polynomial.Monic.irreducible_of_irreducible_map_of_isPrime_nilradical
rw [hm, coeff_mul, Finset.Nat.sum_antidiagonal_eq_sum_range_succ fun i j ↦ a.coeff i * b.coeff j,
Finset.sum_range_succ, ← sub_eq_iff_eq_add, Nat.sub_self] at hc
rw [← add_sub_cancel' 1 (-(_ * _)), ← sub_eq_add_neg, hc]
exact IsNilpotent.sub_one_isUnit <| show _ ∈ nilradical R from sum_mem fun i hi ↦
exact IsNilpotent.isUnit_sub_one <| show _ ∈ nilradical R from sum_mem fun i hi ↦
Ideal.mul_mem_left _ _ <| hn _ <| Nat.sub_ne_zero_of_lt (List.mem_range.1 hi)
31 changes: 26 additions & 5 deletions Mathlib/RingTheory/Polynomial/Nilpotent.lean
Expand Up @@ -6,7 +6,9 @@ Authors: Emilie Uthaiwat, Oliver Nash
import Mathlib.RingTheory.Nilpotent
import Mathlib.Data.Polynomial.AlgebraMap
import Mathlib.Data.Polynomial.Div
import Mathlib.Data.Polynomial.Identities
import Mathlib.RingTheory.Ideal.QuotientOperations
import Mathlib.RingTheory.Polynomial.Tower

/-!
# Nilpotency in polynomial rings.
Expand Down Expand Up @@ -109,11 +111,9 @@ theorem isUnit_of_coeff_isUnit_isNilpotent (hunit : IsUnit (P.coeff 0))
exact hunit.map C }
set P₁ := P.eraseLead with hP₁
suffices IsUnit P₁ by
rw [← eraseLead_add_monomial_natDegree_leadingCoeff P, ← C_mul_X_pow_eq_monomial]
obtain ⟨Q, hQ⟩ := this
rw [← hP₁, ← hQ]
refine' Commute.IsNilpotent.add_isUnit (isNilpotent_C_mul_pow_X_of_isNilpotent _ (hnil _ hdeg))
((Commute.all _ _).mul_left (Commute.all _ _))
rw [← eraseLead_add_monomial_natDegree_leadingCoeff P, ← C_mul_X_pow_eq_monomial, ← hP₁]
refine IsNilpotent.isUnit_add_left_of_commute ?_ this (Commute.all _ _)
exact isNilpotent_C_mul_pow_X_of_isNilpotent _ (hnil _ hdeg)
have hdeg₂ := lt_of_le_of_lt P.eraseLead_natDegree_le (Nat.sub_lt
(Nat.pos_of_ne_zero hdeg) zero_lt_one)
refine' hind P₁.natDegree _ _ (fun i hi => _) rfl
Expand Down Expand Up @@ -171,4 +171,25 @@ lemma isUnit_iff' :

end CommRing

section CommAlgebra

variable {R S : Type*} [CommRing R] [CommRing S] [Algebra R S] (P : R[X]) {a b : S}

lemma isNilpotent_aeval_sub_of_isNilpotent_sub (h : IsNilpotent (a - b)) :
IsNilpotent (aeval a P - aeval b P) := by
simp only [← eval_map_algebraMap]
have ⟨c, hc⟩ := evalSubFactor (map (algebraMap R S) P) a b
exact hc ▸ (Commute.all _ _).isNilpotent_mul_right h

variable {P}

lemma isUnit_aeval_of_isUnit_aeval_of_isNilpotent_sub
(hb : IsUnit (aeval b P)) (hab : IsNilpotent (a - b)) :
IsUnit (aeval a P) := by
rw [← add_sub_cancel'_right (aeval b P) (aeval a P)]
refine IsNilpotent.isUnit_add_left_of_commute ?_ hb (Commute.all _ _)
exact isNilpotent_aeval_sub_of_isNilpotent_sub P hab

end CommAlgebra

end Polynomial
5 changes: 5 additions & 0 deletions Mathlib/RingTheory/Polynomial/Tower.lean
Expand Up @@ -41,6 +41,11 @@ theorem aeval_map_algebraMap (x : B) (p : R[X]) : aeval x (map (algebraMap R A)
rw [aeval_def, aeval_def, eval₂_map, IsScalarTower.algebraMap_eq R A B]
#align polynomial.aeval_map_algebra_map Polynomial.aeval_map_algebraMap

@[simp]
lemma eval_map_algebraMap (P : R[X]) (a : A) :
(map (algebraMap R A) P).eval a = aeval a P := by
rw [← aeval_map_algebraMap (A := A), coe_aeval_eq_eval]

end Semiring

section CommSemiring
Expand Down

0 comments on commit 229b0ba

Please sign in to comment.