Quartic equation

Stream: mathlib4

[Violeta Hernández]

Do we have any progress towards Freek's 46th problem? i.e. the general solution to a quartic

[Violeta Hernández]

It feels by far like the easiest remaining result in that list. Surely you just write everything down and let linarith and ring do the rest?

[Edward van de Meent]

in what form would we want to formalize this? at first glance, i see these two approaches...

import Mathlib.Algebra.Polynomial.Degree.Definitions

import Mathlib.Data.Complex.Basic

set_option autoImplicit false

example (x:ℂ) (a b c d e : ℂ): a * x ^ 4 + b * x ^3 + c * x ^ 2 + d * x + e = 0 → x = sorry :=
  sorry

example (x:ℂ) (p : Polynomial ℂ) (hp : p.degree ≤ 4): p.aeval x = 0 → x = sorry :=
  sorry

[Violeta Hernández]

Looking through other proof assistants' formalizations, it seems they take the first approach.

[Violeta Hernández]

Here's the Coq statement for instance

Theorem Ferrari_formula: forall (a:C) (b:C) (c:C) (d:C),
  let s:=-a/4 in
  let p:= b+3*s*a+6*Cpow s 2 in
  let q:= c+2*b*s+3*a*Cpow s 2+4*Cpow s 3 in
  let r:= d+c*s+b*Cpow s 2+a*Cpow s 3+Cpow s 4 in
  let lambda:=Cardan_Tartaglia_formula (-p/2) (-r) (r*p/2-/8*Cpow q 2) 0 in
  let A:=Csqrt(2*lambda-p) in
  let cond:=(Ceq_dec (2*lambda) p) in
  let B:=if cond then (RtoC 0) else (-q/(2*A)) in
  let z1:=if cond then Csqrt (binom_solution p r 0)
                  else binom_solution A (B+lambda) 0 in
  let z2:=if cond then -Csqrt (binom_solution p r 0)
                  else binom_solution A (B+lambda) 1 in
  let z3:=if cond then Csqrt (binom_solution p r 1)
                  else binom_solution (-A) (-B+lambda) 0 in
  let z4:=if cond then -Csqrt (binom_solution p r 1)
                  else binom_solution (-A) (-B+lambda) 1 in
  let u1:=z1+s in
  let u2:=z2+s in
  let u3:=z3+s in
  let u4:=z4+s in
  forall u:C, (u-u1)*(u-u2)*(u-u3)*(u-u4)=Cpow u 4+a*Cpow u 3+b*Cpow u 2+c*u+d.

[Edward van de Meent]

i'm guessing if we'd ever want to do this, we would like to pull out those lets?

[Violeta Hernández]

Yep, all of those can be auxiliary definitions within the file.

[Violeta Hernández]

By the way, have we done this for the cubic equation already? If so we can just copy their design over

[Ruben Van de Velde]

I think there's three implementations of the cubic in mathlib

[Violeta Hernández]

Can you link to them?

[Edward van de Meent]

i couldn't find a cubic equation solution formula with leansearch...

[Ruben Van de Velde]

Not easily, on mobile. I think one was in archive

[Violeta Hernández]

Oh here it is

[Violeta Hernández]

https://github.com/leanprover-community/mathlib4/blob/ba3db161926d69d9d7b6629b38a98630a32204bb/Archive/Wiedijk100Theorems/SolutionOfCubic.lean#L88-L111

[Violeta Hernández]

I'm not sure what I was expecting but it definitely didn't look like that

[Violeta Hernández]

I guess it all boils down to generality, you don't need to assume the existence of general square and cube roots within your field, only those of certain quantities.

[Violeta Hernández]

...I think the hardest part about formalizing the quartic will be even finding the correct statement

[Edward van de Meent]

i think stating it for complex numbers should do for most applications? (i.e. for applications where there is a ring-embedding into the complex numbers or something like that?)

[Violeta Hernández]

I feel like it'd be a bit of a letdown to have this very general statement for a cubic and just prove the quartic for the specific case of complex numbers

[Violeta Hernández]

I'm sure if we actually dig our way through the formula we can find some similar formulation to the one for the cubic

[Violeta Hernández]

I don't currently have that sort of free time, but maybe in a week or so?

[Edward van de Meent]

in a sense the statement for the cubic maybe isn't that general either? you don't need every number to have an inverse afaict...

[Violeta Hernández]

Cubics over a general ring need not have only three solutions!

[Edward van de Meent]

over domains they do

[Edward van de Meent]

see also Polynomial.card_roots

[Violeta Hernández]

I wonder what breaks if you try proving the cubic equation theorem for a commutative domain

[Violeta Hernández]

Like, surely you can still substitute in the values and have them work? And then argue that there can be no other roots

[Edward van de Meent]

the values don't need to exist in the domain, and i don't know that you can then still rule out the existance of replacement "unnatural" solutions?

[Violeta Hernández]

What do you mean? The values are defined in terms of sums, products, and divisions by invertible constants

[Edward van de Meent]

aren't there any cube roots?

[Violeta Hernández]

nope!

[Violeta Hernández]

that's the nice thing about this formulation

[Edward van de Meent]

still, the constants needn't be invertible?

[Violeta Hernández]

theorem cubic_eq_zero_iff (ha : a ≠ 0) (hω : IsPrimitiveRoot ω 3)
    (hp : p = (3 * a * c - b ^ 2) / (9 * a ^ 2)) (hp_nonzero : p ≠ 0)
    (hq : q = (9 * a * b * c - 2 * b ^ 3 - 27 * a ^ 2 * d) / (54 * a ^ 3))
    (hr : r ^ 2 = q ^ 2 + p ^ 3) (hs3 : s ^ 3 = q + r) (ht : t * s = p) (x : K) :
    a * x ^ 3 + b * x ^ 2 + c * x + d = 0 ↔
      x = s - t - b / (3 * a) ∨
        x = s * ω - t * ω ^ 2 - b / (3 * a) ∨ x = s * ω ^ 2 - t * ω - b / (3 * a) := sorry

[Violeta Hernández]

Yeah, there's a typeclass assumption that 2 and 3 are invertible

[Violeta Hernández]

But that's all you really need

[Violeta Hernández]

If there's a statement like this for the quartic equation then that's probably the cleanest way to prove it

[Violeta Hernández]

Oh wait of course you need a to be invertible too

[Violeta Hernández]

that's just being covered by a ≠ 0 in this case

[Violeta Hernández]

The other nice thing about this formulation is that you don't have to argue about the multivalued nature of square and cube roots. In the textbook formulation of the cubic formula, you have this weirdness where there's 6 possible values but only three of them are distinct. Introducing auxiliary variables like this gets rid of this issue

[Edward van de Meent]

makes sense...

[Violeta Hernández]

I'll definitely keep this project in mind. No promises though :stuck_out_tongue:

[Edward van de Meent]

i'm giving it a small try right now :upside_down:

[Kevin Buzzard]

Every integral domain is a subring of a field so solving in integral domains is basically the same as solving in fields

[Thomas Zhu]

While searching on google I found this:
https://github.com/anonymousLeanDocsHosting/lean-polynomials/blob/main/mathlib/quartic_solution.lean
and its paper. Apparently it was already done in Lean 3, assuming it builds correctly?

[Thomas Zhu]

By the way, the Isabelle/HOL proof of the quartic formula is very short (4 lines), though it might have assumed that coefficients are real? It's the same formulation as the one listed on Wolfram

[Thomas Zhu]

I just had a go at the quartic formula, and I have proven it:

namespace Theorems100

open Polynomial

section Quartic

variable {K : Type*} [Field K] (a b c d e : K) {ω p q r s t u v w x y z : K}

variable [Invertible (2 : K)]

/-- Roots of a quartic whose cubic term is zero and linear term is nonzero. -/
theorem quartic_basic_eq_zero_iff
    (hq_nonzero : q ≠ 0)
    -- u satisfies the cubic resolvent
    (hu : u ^ 3 - p * u ^ 2 - 4 * r * u + 4 * p * r - q ^ 2 = 0)
    (hs : s ^ 2 = u - p)
    (hv : v ^ 2 = 4 * s ^ 2 - 8 * (u - q / s))
    (hw : w ^ 2 = 4 * s ^ 2 - 8 * (u + q / s))
    (x : K) :
    x ^ 4 + p * x ^ 2 + q * x + r = 0 ↔
      x = (-2 * s - v) / 4 ∨ x = (-2 * s + v) / 4 ∨ x = (2 * s - w) / 4 ∨ x = (2 * s + w) / 4 := by
  have hi2 : (2 : K) ≠ 0 := Invertible.ne_zero _
  have h4 : (4 : K) = 2 ^ 2 := by norm_num
  have hs_nonzero : s ≠ 0 := by
    contrapose! hq_nonzero with hs0
    linear_combination (exp := 2) -hu + (4 * r - u ^ 2) * hs + (u ^ 2 * s - 4 * r * s) * hs0
  calc
    _ ↔ 4 * (x ^ 4 + p * x ^ 2 + q * x + r) = 0 := by simp [h4, hi2]
    _ ↔ (2 * (x * x) + 2 * s * x + (u - q / s)) * (2 * (x * x) + -(2 * s) * x + (u + q / s)) =
        0 := by
      refine Eq.congr ?_ rfl
      field_simp
      linear_combination -hu + (-x ^ 2 * s ^ 2 - x ^ 2 * p + x ^ 2 * u) * hw +
        (x ^ 2 * w ^ 2 + 8 * x ^ 2 * u + 8 * x ^ 2 * q / s - u ^ 2 + 4 * r) * hs
    _ ↔ _ := by
      have hv' : discrim 2 (2 * s) (u - q / s) = v * v := by rw [discrim]; linear_combination -hv
      have hw' : discrim 2 (-(2 * s)) (u + q / s) = w * w := by rw [discrim]; linear_combination -hw
      rw [mul_eq_zero, quadratic_eq_zero_iff hi2 hv', quadratic_eq_zero_iff hi2 hw']
      simp [(by norm_num : (2 : K) * 2 = 4), or_assoc, or_comm]

/-- **The Solution of Quartic**.
  The roots of a quartic polynomial when q is nonzero.
  See Beyer W. H., CRC Standard Mathematical Tables and Formulae.
  Here, u needs to satisfy the cubic resolvent. An explicit expression of u is possible using
  the cubic formula, but would be too long. -/
theorem quartic_eq_zero_iff (ha : a ≠ 0)
    (hp : p = (8 * a * c - 3 * b ^ 2) / (8 * a ^ 2))
    (hq : q = (b ^ 3 - 4 * a * b * c + 8 * a ^ 2 * d) / (8 * a ^ 3)) (hq_nonzero : q ≠ 0)
    (hr : r =
      (16 * a * b ^ 2 * c + 256 * a ^ 3 * e - 3 * b ^ 4 - 64 * a ^ 2 * b * d) / (256 * a ^ 4))
    (hu : u ^ 3 - p * u ^ 2 - 4 * r * u + 4 * p * r - q ^ 2 = 0)
    (hs : s ^ 2 = u - p)
    (hv : v ^ 2 = 4 * s ^ 2 - 8 * (u - q / s))
    (hw : w ^ 2 = 4 * s ^ 2 - 8 * (u + q / s)) (x : K) :
    a * x ^ 4 + b * x ^ 3 + c * x ^ 2 + d * x + e = 0 ↔
      x = (-2 * s - v) / 4 - b / (4 * a) ∨ x = (-2 * s + v) / 4 - b / (4 * a) ∨
        x = (2 * s - w) / 4 - b / (4 * a) ∨ x = (2 * s + w) / 4 - b / (4 * a) := by
  let y := x + b / (4 * a)
  have h4 : (4 : K) = 2 ^ 2 := by norm_num
  have h8 : (8 : K) = 2 ^ 3 := by norm_num
  have h16 : (16 : K) = 2 ^ 4 := by norm_num
  have h256 : (256 : K) = 2 ^ 8 := by norm_num
  have h₁ : a * x ^ 4 + b * x ^ 3 + c * x ^ 2 + d * x + e =
      a * (y ^ 4 + p * y ^ 2 + q * y + r) := by
    rw [hp, hq, hr]
    simp [field_simps, y, ha, h4, h8, h16, h256]; ring
  have h₂ : ∀ x, a * x = 0 ↔ x = 0 := by intro x; simp [ha]
  rw [h₁, h₂, quartic_basic_eq_zero_iff hq_nonzero hu hs hv hw]
  simp_rw [eq_sub_iff_add_eq]

/-- The roots of a quartic polynomial when q equals zero. -/
theorem quartic_eq_zero_iff_of_zero (ha : a ≠ 0)
    (hp : p = (8 * a * c - 3 * b ^ 2) / (8 * a ^ 2))
    (hqz : b ^ 3 - 4 * a * b * c + 8 * a ^ 2 * d = 0)
    (hr : r =
      (16 * a * b ^ 2 * c + 256 * a ^ 3 * e - 3 * b ^ 4 - 64 * a ^ 2 * b * d) / (256 * a ^ 4))
    (ht : t ^ 2 = p ^ 2 - 4 * r)
    (hv : v ^ 2 = (-p + t) / 2)
    (hw : w ^ 2 = (-p - t) / 2) (x : K) :
    a * x ^ 4 + b * x ^ 3 + c * x ^ 2 + d * x + e = 0 ↔
      x = v - b / (4 * a) ∨ x = -v - b / (4 * a) ∨ x = w - b / (4 * a) ∨ x = -w - b / (4 * a) := by
  let y := x + b / (4 * a)
  have h4 : (4 : K) = 2 ^ 2 := by norm_num
  have h8 : (8 : K) = 2 ^ 3 := by norm_num
  have h16 : (16 : K) = 2 ^ 4 := by norm_num
  have h256 : (256 : K) = 2 ^ 8 := by norm_num
  have h₁ : a * x ^ 4 + b * x ^ 3 + c * x ^ 2 + d * x + e = a * (y ^ 4 + p * y ^ 2 + r) := by
    rw [hp, hr]
    simp [field_simps, y, ha, h4, h8, h16, h256]
    linear_combination (1048576 * a ^ 10 * x + 262144 * a ^ 9 * b) * hqz
  have h₂ : ∀ x, a * x = 0 ↔ x = 0 := by intro x; simp [ha]
  rw [h₁, h₂]
  calc
    _ ↔ 1 * (y ^ 2 * y ^ 2) + p * y ^ 2 + r = 0 := by
      refine Eq.congr ?_ rfl
      ring
    _ ↔ y ^ 2 = (-p + t) / 2 ∨ y ^ 2 = (-p - t) / 2 := by
      have ht' : discrim 1 p r = t * t := by rw [discrim]; linear_combination -ht
      rw [quadratic_eq_zero_iff one_ne_zero ht', mul_one]
    _ ↔ _ := by
      simp_rw [← hv, ← hw, pow_two, mul_self_eq_mul_self_iff, eq_sub_iff_add_eq, or_assoc]

end Quartic

end Theorems100

[Thomas Zhu]

The basic idea is the same as cubic: you can reduce a general quartic polynomial $ax^4+bx^3+cx^2+dx+e=0$ to one with cubic coefficient = 0 (the depressed/basic form), here $y ^ 4 + p y ^ 2 + q y + r=0$. In this form, either the $q$ is nonzero, in which case there is a further reduction to solving a cubic polynomial (the cubic resolvent), or $q=0$ in which case it is just a quadratic in $y^2$. The formulation of the cubic resolvent differs in different statements of the theorem; I have chosen a particular one.

[Thomas Zhu]

PR in progress! #18290