Defining n-division polynomials

Stream: maths

[Marc Masdeu]

Hi all,

I'm having a go at proving in Lean something that can't be found in the literature (as far as I'm aware), namely a purely algebraic proof of the validity of the formulas for the multiplication-by-n isogeny on elliptic curves. This involves the infamous n-division polynomials, which I am trying to define in Lean. The recurrence has 5 basic cases (0,1,2,3,4) and then two recurrent formulas, depending on the parity of n. I haven't been able to make it work so far, here is my try:

import tactic
import data.polynomial.basic
open polynomial

noncomputable theory
--variables {R : Type} [comm_ring R]
variables {A B : ℤ}
def x := (monomial 1 (1 : ℤ))

--set_option trace.eqn_compiler.elim_match true
noncomputable def psi : ℕ → (polynomial ℤ)
| 0 := 0
| 1 := 1
| 2 := 1
| 3 := 3*x^4 + 6*A*x^2 + 12*B*x - A^2*x
| 4 := 2*(x^6 + 5*A*x^4+20*B*x^3-5*A^2*x^2-4*A*B*x-A^3-8*B^2)
| (2*n+5) :=
  have n+4 < 2*n+5, by linarith,
  have n+3 < 2*n+5, by linarith,
  have n+2 < 2*n+5, by linarith,
  have n+1 < 2*n+5, by linarith,
  (psi (n+4)) * (psi (n+2))^3 - (psi (n+1) * (psi (n+3))^3)
| (2*n+6) :=
  have n+5 < 2*n+6, by linarith,
  have n+4 < 2*n+6, by linarith,
  have n+3 < 2*n+6, by linarith,
  have n+2 < 2*n+6, by linarith,
  (psi (n+3)) * ((psi (n+5))*
      (psi (n+2))^2 - (psi (n+1)) * (psi (n+4))^2) / (4*(x^3+A*x+B)),

From reading https://leanprover-community.github.io/extras/well_founded_recursion.html I gathered that I needed to supply the proofs of the fact that the recurrence is well-founded, but they would all be easy, I could even figure out how to produce a term proof if linarith was forbidden.

Also, I'd like A and B to be allowed to live in any commutative ring, but then it gives me even more errors because apparently A*x is not a polynomial in R when A is in R (but that's for a different thread I guess).

Any hints would be appreciated :-).

PS: feel free to move this to "new members", I'm not always sure of the distinction between streams.

[Kevin Buzzard]

I think you're no longer new :-) Just to be clear -- this is the y^2=x^3+Ax+B model, right?

[Kevin Buzzard]

The equation compiler isn't going to be able to split between even and odd cases I don't think, because even v odd is not something which comes directly from the recursive definition of the naturals. I would just be tempted to use "strong recursion" here.

[Marc Masdeu]

(what if not new but have a newbie question?). This is the y^2=x^3+Ax+B model, but if that worked then it shouldn't be too hard to adapt the proof to the general Weierstrass model.

[Marc Masdeu]

Oh, so you say that 2n vs 2n+1 (or what I write) does confuse lean? Maybe there's an easier way just defining phi_odd and phi_even, then? It wouldn't be as elegant, but I just want to type a proof...

[Johan Commelin]

Marc Masdeu said:

Also, I'd like A and B to be allowed to live in any commutative ring, but then it gives me even more errors because apparently A*x is not a polynomial in R when A is in R (but that's for a different thread I guess).

Use C A * X or A \bullet X. Here C is polynomial.C is the ring morphism from R to polynomial R.

(I think this thread perfectly fits the #maths stream!)

[Johan Commelin]

@Marc Masdeu I think would use if .. then .. else in the n + 5 branch

[Johan Commelin]

And nat-division, which is truncating, so 5/2 = 2, 6/2 = 3, and 7/2 = 3, etc...

[Johan Commelin]

The other option is the strong recursion that Kevin mentioned

[Kevin Buzzard]

Aah very nice, so you don't have to extract the k from "exists k, n=2k or n=2k+1"

[Kevin Buzzard]

I was going to suggest matching on the pos_num :-) but integer division is a much more sensible approach

[Kevin Buzzard]

So in conclusion do the final match on n+5 and use if then else on even (n+5) (or probably even n is more sensible) and try your have technique. And for the polynomials I would be tempted to fix a ring and just use polynomial.X directly (open polynomial first) and then use C or bub

[Marc Masdeu]

Thanks for the suggestions, I'll try and see.

[Marc Masdeu]

Seems that I've made a bit of progress:

import tactic
import data.polynomial.basic
open polynomial

noncomputable theory
variables {R : Type} [comm_ring R]
variables {A B : R}
def x := (monomial 1 (1 : R))

--set_option trace.eqn_compiler.elim_match true
noncomputable def psi : ℕ → (polynomial R)
| 0 := 0
| 1 := 1
| 2 := 1
| 3 := 3*x^4 + 6*A•x^2 + 12*B•x - A^2•x
| 4 := 2*(x^6 + 5*A•x^4+20*B•x^3-5*A^2•x^2-(4*A*B)•x-A^3•x^0-8*B^2•x^0)
| (n+5) :=
  have n/2+1 < n+5, by sorry,
  have n/2+2 < n+5, by sorry,
  have n/2+3 < n+5, by sorry,
  have n/2+4 < n+5, by sorry,
  have n/2+5 < n+5, by sorry,
  if even n then (psi (n/2+3)) * ((psi (n/2+5))*
      (psi (n/2+2))^2 - (psi (n/2+1)) * (psi (n/2+4))^2) --/ (4*(x^3+A•x+B•x^0))
      else
  (psi (n/2+4)) * (psi (n/2+2))^3 - (psi (n/2+1) * (psi (n/2+3))^3)

but Lean complains about decidable (even n). Also I'll have trouble dividing by a polynomial, I need to familiarize myself with that part of mathlib, but I'd have assumed that there is a has_div giving out junk when the numerator is not a multiple of the denominator, as usual.

[Johan Commelin]

I'm surprised that it complains about decidable (even n). But I guess you can replace it with n % 2 = 0 or open_locale classical.

[Johan Commelin]

n % 2 = 0 might be a bit easier, because it allows you to use rfl in a bit more places.

[Johan Commelin]

Of course the other option is to write an instance for decidable (even n). I didn't find it in mathlib with a quick search.

[Johan Commelin]

Anyway, the first thing after this definition is probably simp-lemmas for 0, 1, 2, 3, 4, even n, and odd n. Maybe also for 2n + 5 and 2n + 6.

[Kevin Buzzard]

And also I would be tempted to use polynomial.X rather than your x because there will be a bunch of lemmas about how X evaluates but you'll be constantly unfolding your x in order to make progress

[Kevin Buzzard]

import tactic
import data.nat.parity

example (n : ℕ) : decidable (even n) := infer_instance

It is there, you just need the right import perhaps?

[Marc Masdeu]

For what is worth, I post the definition I settled on.

import tactic
import data.polynomial.basic
import data.polynomial.div
import data.nat.parity
open polynomial
open nat

noncomputable theory
open_locale classical

lemma nhalves {n m r : ℕ} :  n / 2 + m < n + m + r + 1 :=
begin
  have twopos : 0 < 2, by linarith,
  have : n < 2*(n+1), by linarith,
  suffices : n / 2  < n+1,
  by linarith,
  exact nat.div_lt_of_lt_mul this,
end

lemma nhalves' {n m r : ℕ} : (n+1)/2 + m < n + m + r + 1:=
begin
  suffices : (n+1)/2 + m < n + m + 1, by linarith,
  suffices : (n+1)/2 < n + 1, by linarith,
  exact nat.div_lt_self' n 0,
end

/--
Define the "depleted" psi polynomials. They are Ψ as in Wikipedia,
but for even n we remove a factor of 2*y, thus making them one-variable.
-/
noncomputable def psi {R : Type*} [comm_ring R] (A B : R) : ℕ → (polynomial R)
| 0 := (0 : polynomial R)
| 1 := (1 : polynomial R)
| 2 := (1 : polynomial R)
| 3 := 3*X^4 + 6*A•X^2 + 12*B•X - A^2•1
| 4 := 2*(X^6 + 5*A•X^4+20*B•X^3-5*A^2•X^2-(4*A*B)•X-A^3•1-8*B^2•1)
| (n+5) :=
  have n/2+4 < n+4+0+1 := nhalves,
  have n/2+3 < n+3+1+1 := nhalves,
  have n/2+2 < n+2+2+1 := nhalves,
  have n/2+1 < n+1+3+1 := nhalves,
  have (n+1)/2 < n+0+4+1 := @nhalves' _ 0 _,
  have (n+1)/2+1 < n+1+3+1 := nhalves',
  have (n+1)/2+2 < n+2+2+1 := nhalves',
  have (n+1)/2+3 < n+3+1+1 := nhalves',
  have (n+1)/2+4 < n+4+0+1 := nhalves',
  if even n then
    if even (n / 2) then
      16*(X^3+A•X+B•1)^2 * psi (n/2+4) * (psi (n/2+2))^3 - psi (n/2+1) * (psi (n/2+3))^3
    else
      psi (n/2+4) * (psi (n/2+2))^3 - 16*(X^3+A•X+B•1)^2 * psi (n/2+1) * (psi (n/2+3))^3
  else
      (psi ((n+1)/2+2)) * ((psi ((n+1)/2+4)) * (psi ((n+1)/2+1))^2 -
        psi ((n+1)/2) * (psi ((n+1)/2+3))^2)

I have been able to prove the lemma, for example, with quite a bit of work. But it would be useful to avoid having to write A B so many times in such statements. For my application, I'd like to "fix" A B throughout a namespace, and define something like mypsi (or better, Ψ) like Ψ : ℕ → polynomial R := λ n, psi A B n, but I don't know how to do it, really.

lemma psi_two_times {A B : R} {n : ℕ}  (h : 2 ≤ n) :
 psi A B (2*n) = (psi A B n) * (psi A B (n+2) * (psi A B (n-1))^2 -
  psi A B (n-2) * (psi A B (n+1))^2 ) := sorry

[Johan Commelin]

@Marc Masdeu One possibility would be to use notation for this:

local notation `mypsi` := `psi A B`

That way you don't have to write the A and B everywhere in your code. But they will still be visible in the goal window.

[Mario Carneiro]

you mean local notation `mypsi` := psi A B

[Damiano Testa]

Hi Marc!

Not at all related to the A B stuff, but the initial lemmas can be proved as follows. I was actually surprised that lemma nat.half_le_self is not already in mathlib (and maybe it is?):

lemma nat.half_le_self (n : ℕ) : (n + 1) / 2 ≤ n :=
begin
  cases n,
  { exact rfl.le },
  { exact nat.div_le_of_le_mul (((add_le_add_iff_left _).mpr ((1 : ℕ).le_add_left n)).trans
      (two_mul _).symm.le) }
end

lemma nhalves {n m r : ℕ} :  n / 2 + m < n + m + r + 1 :=
(((add_le_add_iff_right _).mpr (n.div_le_self 2)).trans (le.intro rfl)).trans_lt (lt_add_one _)

lemma nhalves' {n m r : ℕ} : (n+1)/2 + m < n + m + r + 1:=
(((add_le_add_iff_right _).mpr n.half_le_self).trans (le.intro rfl)).trans_lt (lt_add_one _)

[Yakov Pechersky]

Check if there is such a statement about div2

[Damiano Testa]

I thought that I checked: besides library_search, I also tried with the autocomplete, but all my guesses failed. I did find the one with n / 2 ≤ n, though...

[Damiano Testa]

There is a possibly promising lemma, but I am not sure what it says:

@[simp] lemma div2_succ (n : ℕ) : div2 (succ n) = cond (bodd n) (succ (div2 n)) (div2 n) :=
by unfold bodd div2 bodd_div2; cases bodd_div2 n; cases fst; refl

[Damiano Testa]

Also, I am at the same time proud and ashamed of having guessed that (0 + 1) / 2 ≤ 0 was

• true,
• proved by exact rfl.le.

[Damiano Testa]

Johan, I wish I could add emoji reactions to your emoji reactions! :lol:

[Damiano Testa]

In case you find this helpful, I made a PR with the lemma (n+1) / 2 ≤ n.
#6863

[Damiano Testa]

Marc, lemma nat.succ_half_le_self did not survive the review process. However, with the latest version of mathlib, what's below is also a proof.

lemma nat.succ_half_le_self (n : ℕ) : (n + 1) / 2 ≤ n :=
(div_le_iff_le_mul_add_pred zero_lt_two).mpr ((add_le_add_iff_right _).mpr
  (le_mul_of_one_le_of_le one_le_two rfl.le))

If you do not care about lines longer than 100 chars, it even fits on a single line. :upside_down:

[Bhavik Mehta]

Every proof is a one line proof if you start far enough to the left :upside_down:

[Johan Commelin]

By induction every proof fits in the margin?