(Joint Center)Library MathComp.prime

(* (c) Copyright Microsoft Corporation and Inria.                       
 You may distribute this file under the terms of the CeCILL-B license *)

Require Import ssreflect ssrfun ssrbool eqtype ssrnat seq path fintype.
Require Import div bigop.

This file contains the definitions of: prime p <=> p is a prime. primes m == the sorted list of prime divisors of m > 1, else [:: ]. pfactor == the type of prime factors, syntax (p ^ e)%pfactor. prime_decomp m == the list of prime factors of m > 1, sorted by primes. logn p m == the e such that (p ^ e) \in prime_decomp n, else 0. trunc_log p m == the largest e such that p ^ e <= m, or 0 if p or m is 0. pdiv n == the smallest prime divisor of n > 1, else 1. max_pdiv n == the largest prime divisor of n > 1, else 1. divisors m == the sorted list of divisors of m > 0, else [:: ]. totient n == the Euler totient (#|{i < n | i and n coprime}|). nat_pred == the type of explicit collective nat predicates. := simpl_pred nat.
  • > We allow the coercion nat >-> nat_pred, interpreting p as pred1 p.
  • > We define a predType for nat_pred, enabling the notation p \in pi.
  • > We don't have nat_pred >-> pred, which would imply nat >-> Funclass. pi^' == the complement of pi : nat_pred, i.e., the nat_pred such that (p \in pi^') = (p \notin pi). \pi(n) == the set of prime divisors of n, i.e., the nat_pred such that (p \in \pi(n)) = (p \in primes n). \pi(A) == the set of primes of #|A|, with A a collective predicate over a finite Type.
    • > The notation \pi(A) is implemented with a collapsible Coercion, so the type of A must coerce to finpred_class (e.g., by coercing to {set T}), not merely implement the predType interface (as seq T does).
    • > The expression #|A| will only appear in \pi(A) after simplification collapses the coercion stack, so it is advisable to do so early on.
    pi.-nat n <=> n > 0 and all prime divisors of n are in pi. n`_pi == the pi-part of n -- the largest pi.-nat divisor of n. := \prod_(0 <= p < n.+1 | p \in pi) p ^ logn p n.
    • > The nat >-> nat_pred coercion lets us write p.-nat n and n`_p.
In addition to the lemmas relevant to these definitions, this file also contains the dvdn_sum lemma, so that bigop.v doesn't depend on div.v.

Set Implicit Arguments.

The complexity of any arithmetic operation with the Peano representation is pretty dreadful, so using algorithms for "harder" problems such as factoring, that are geared for efficient artihmetic leads to dismal performance -- it takes a significant time, for instance, to compute the divisors of just a two-digit number. On the other hand, for Peano integers, prime factoring (and testing) is linear-time with a small constant factor -- indeed, the same as converting in and out of a binary representation. This is implemented by the code below, which is then used to give the "standard" definitions of prime, primes, and divisors, which can then be used casually in proofs with moderately-sized numeric values (indeed, the code here performs well for up to 6-digit numbers).
We start with faster mod-2 functions.

Fixpoint edivn2 q r := if r is r'.+2 then edivn2 q.+1 r' else (q, r).

Lemma edivn2P n : edivn_spec n 2 (edivn2 0 n).

Fixpoint elogn2 e q r {struct q} :=
  match q, r with
  | 0, _ | _, 0 ⇒ (e, q)
  | q'.+1, 1 ⇒ elogn2 e.+1 q' q'
  | q'.+1, r'.+2elogn2 e q' r'
  end.

CoInductive elogn2_spec n : nat × natType :=
  Elogn2Spec e m of n = 2 ^ e × m.*2.+1 : elogn2_spec n (e, m).

Lemma elogn2P n : elogn2_spec n.+1 (elogn2 0 n n).

Definition ifnz T n (x y : T) := if n is 0 then y else x.

CoInductive ifnz_spec T n (x y : T) : TType :=
  | IfnzPos of n > 0 : ifnz_spec n x y x
  | IfnzZero of n = 0 : ifnz_spec n x y y.

Lemma ifnzP T n (x y : T) : ifnz_spec n x y (ifnz n x y).

For pretty-printing.
Definition NumFactor (f : nat × nat) := ([Num of f.1], f.2).

Definition pfactor p e := p ^ e.

Definition cons_pfactor (p e : nat) pd := ifnz e ((p, e) :: pd) pd.

Notation Local "p ^? e :: pd" := (cons_pfactor p e pd)
  (at level 30, e at level 30, pd at level 60) : nat_scope.

Section prime_decomp.

Import NatTrec.

Fixpoint prime_decomp_rec m k a b c e :=
  let p := k.*2.+1 in
  if a is a'.+1 then
    if b - (ifnz e 1 k - c) is b'.+1 then
      [rec m, k, a', b', ifnz c c.-1 (ifnz e p.-2 1), e] else
    if (b == 0) && (c == 0) then
      let b' := k + a' in [rec b'.*2.+3, k, a', b', k.-1, e.+1] else
    let bc' := ifnz e (ifnz b (k, 0) (edivn2 0 c)) (b, c) in
    p ^? e :: ifnz a' [rec m, k.+1, a'.-1, bc'.1 + a', bc'.2, 0] [:: (m, 1)]
  else if (b == 0) && (c == 0) then [:: (p, e.+2)] else p ^? e :: [:: (m, 1)]
where "[ 'rec' m , k , a , b , c , e ]" := (prime_decomp_rec m k a b c e).

Definition prime_decomp n :=
  let: (e2, m2) := elogn2 0 n.-1 n.-1 in
  if m2 < 2 then 2 ^? e2 :: 3 ^? m2 :: [::] else
  let: (a, bc) := edivn m2.-2 3 in
  let: (b, c) := edivn (2 - bc) 2 in
  2 ^? e2 :: [rec m2.*2.+1, 1, a, b, c, 0].

The list of divisors and the Euler function are computed directly from the decomposition, using a merge_sort variant sort the divisor list.

Definition add_divisors f divs :=
  let: (p, e) := f in
  let add1 divs' := merge leq (map (NatTrec.mul p) divs') divs in
  iter e add1 divs.

Definition add_totient_factor f m := let: (p, e) := f in p.-1 × p ^ e.-1 × m.

End prime_decomp.

Definition primes n := unzip1 (prime_decomp n).

Definition prime p := if prime_decomp p is [:: (_ , 1)] then true else false.

Definition nat_pred := simpl_pred nat.

Definition pi_unwrapped_arg := nat.
Definition pi_wrapped_arg := wrapped nat.
Coercion unwrap_pi_arg (wa : pi_wrapped_arg) : pi_unwrapped_arg := unwrap wa.
Coercion pi_arg_of_nat (n : nat) := Wrap n : pi_wrapped_arg.
Coercion pi_arg_of_fin_pred T pT (A : @fin_pred_sort T pT) : pi_wrapped_arg :=
  Wrap #|A|.

Definition pi_of (n : pi_unwrapped_arg) : nat_pred := [pred p in primes n].

Notation "\pi ( n )" := (pi_of n)
  (at level 2, format "\pi ( n )") : nat_scope.
Notation "\p 'i' ( A )" := \pi(#|A|)
  (at level 2, format "\p 'i' ( A )") : nat_scope.

Definition pdiv n := head 1 (primes n).

Definition max_pdiv n := last 1 (primes n).

Definition divisors n := foldr add_divisors [:: 1] (prime_decomp n).

Definition totient n := foldr add_totient_factor (n > 0) (prime_decomp n).

Correctness of the decomposition algorithm.

Lemma prime_decomp_correct :
  let pd_val pd := \prod_(f <- pd) pfactor f.1 f.2 in
  let lb_dvd q m := ~~ has [pred d | d %| m] (index_iota 2 q) in
  let pf_ok f := lb_dvd f.1 f.1 && (0 < f.2) in
  let pd_ord q pd := path ltn q (unzip1 pd) in
  let pd_ok q n pd := [/\ n = pd_val pd, all pf_ok pd & pd_ord q pd] in
   n, n > 0 → pd_ok 1 n (prime_decomp n).

Lemma primePn n :
  reflect (n < 2 exists2 d, 1 < d < n & d %| n) (~~ prime n).

Lemma primeP p :
  reflect (p > 1 d, d %| pxpred2 1 p d) (prime p).

Lemma prime_nt_dvdP d p : prime pd != 1 → reflect (d = p) (d %| p).

Implicit Arguments primeP [p].
Implicit Arguments primePn [n].

Lemma prime_gt1 p : prime p → 1 < p.

Lemma prime_gt0 p : prime p → 0 < p.

Hint Resolve prime_gt1 prime_gt0.

Lemma prod_prime_decomp n :
  n > 0 → n = \prod_(f <- prime_decomp n) f.1 ^ f.2.

Lemma even_prime p : prime pp = 2 odd p.

Lemma prime_oddPn p : prime preflect (p = 2) (~~ odd p).

Lemma odd_prime_gt2 p : odd pprime pp > 2.

Lemma mem_prime_decomp n p e :
  (p, e) \in prime_decomp n[/\ prime p, e > 0 & p ^ e %| n].

Lemma prime_coprime p m : prime pcoprime p m = ~~ (p %| m).

Lemma dvdn_prime2 p q : prime pprime q(p %| q) = (p == q).

Lemma Euclid_dvdM m n p : prime p(p %| m × n) = (p %| m) || (p %| n).

Lemma Euclid_dvd1 p : prime p(p %| 1) = false.

Lemma Euclid_dvdX m n p : prime p(p %| m ^ n) = (p %| m) && (n > 0).

Lemma mem_primes p n : (p \in primes n) = [&& prime p, n > 0 & p %| n].

Lemma sorted_primes n : sorted ltn (primes n).

Lemma eq_primes m n : (primes m =i primes n) (primes m = primes n).

Lemma primes_uniq n : uniq (primes n).

The smallest prime divisor

Lemma pi_pdiv n : (pdiv n \in \pi(n)) = (n > 1).

Lemma pdiv_prime n : 1 < nprime (pdiv n).

Lemma pdiv_dvd n : pdiv n %| n.

Lemma pi_max_pdiv n : (max_pdiv n \in \pi(n)) = (n > 1).

Lemma max_pdiv_prime n : n > 1 → prime (max_pdiv n).

Lemma max_pdiv_dvd n : max_pdiv n %| n.

Lemma pdiv_leq n : 0 < npdiv n n.

Lemma max_pdiv_leq n : 0 < nmax_pdiv n n.

Lemma pdiv_gt0 n : 0 < pdiv n.

Lemma max_pdiv_gt0 n : 0 < max_pdiv n.
Hint Resolve pdiv_gt0 max_pdiv_gt0.

Lemma pdiv_min_dvd m d : 1 < dd %| mpdiv m d.

Lemma max_pdiv_max n p : p \in \pi(n)p max_pdiv n.

Lemma ltn_pdiv2_prime n : 0 < nn < pdiv n ^ 2 → prime n.

Lemma primePns n :
  reflect (n < 2 p, [/\ prime p, p ^ 2 n & p %| n]) (~~ prime n).

Implicit Arguments primePns [n].

Lemma pdivP n : n > 1 → {p | prime p & p %| n}.

Lemma primes_mul m n p : m > 0 → n > 0 →
  (p \in primes (m × n)) = (p \in primes m) || (p \in primes n).

Lemma primes_exp m n : n > 0 → primes (m ^ n) = primes m.

Lemma primes_prime p : prime pprimes p = [::p].

Lemma coprime_has_primes m n : m > 0 → n > 0 →
  coprime m n = ~~ has (mem (primes m)) (primes n).

Lemma pdiv_id p : prime ppdiv p = p.

Lemma pdiv_pfactor p k : prime ppdiv (p ^ k.+1) = p.

"prime" logarithms and p-parts.

Fixpoint logn_rec d m r :=
  match r, edivn m d with
  | r'.+1, (_.+1 as m', 0)(logn_rec d m' r').+1
  | _, _ ⇒ 0
  end.

Definition logn p m := if prime p then logn_rec p m m else 0.

Lemma lognE p m :
  logn p m = if [&& prime p, 0 < m & p %| m] then (logn p (m %/ p)).+1 else 0.

Lemma logn_gt0 p n : (0 < logn p n) = (p \in primes n).

Lemma ltn_log0 p n : n < plogn p n = 0.

Lemma logn0 p : logn p 0 = 0.

Lemma logn1 p : logn p 1 = 0.

Lemma pfactor_gt0 p n : 0 < p ^ logn p n.
Hint Resolve pfactor_gt0.

Lemma pfactor_dvdn p n m : prime pm > 0 → (p ^ n %| m) = (n logn p m).

Lemma pfactor_dvdnn p n : p ^ logn p n %| n.

Lemma logn_prime p q : prime qlogn p q = (p == q).

Lemma pfactor_coprime p n :
  prime pn > 0 → {m | coprime p m & n = m × p ^ logn p n}.

Lemma pfactorK p n : prime plogn p (p ^ n) = n.

Lemma pfactorKpdiv p n : prime plogn (pdiv (p ^ n)) (p ^ n) = n.

Lemma dvdn_leq_log p m n : 0 < nm %| nlogn p m logn p n.

Lemma ltn_logl p n : 0 < nlogn p n < n.

Lemma logn_Gauss p m n : coprime p mlogn p (m × n) = logn p n.

Lemma lognM p m n : 0 < m → 0 < nlogn p (m × n) = logn p m + logn p n.

Lemma lognX p m n : logn p (m ^ n) = n × logn p m.

Lemma logn_div p m n : m %| nlogn p (n %/ m) = logn p n - logn p m.

Lemma dvdn_pfactor p d n : prime p
  reflect (exists2 m, m n & d = p ^ m) (d %| p ^ n).

Lemma prime_decompE n : prime_decomp n = [seq (p, logn p n) | p <- primes n].

Some combinatorial formulae.

Lemma divn_count_dvd d n : n %/ d = \sum_(1 i < n.+1) (d %| i).

Lemma logn_count_dvd p n : prime plogn p n = \sum_(1 k < n) (p ^ k %| n).

Truncated real log.

Definition trunc_log p n :=
  let fix loop n k :=
    if k is k'.+1 then if p n then (loop (n %/ p) k').+1 else 0 else 0
  in loop n n.

Lemma trunc_log_bounds p n :
  1 < p → 0 < nlet k := trunc_log p n in p ^ k n < p ^ k.+1.

Lemma trunc_log_ltn p n : 1 < pn < p ^ (trunc_log p n).+1.

Lemma trunc_logP p n : 1 < p → 0 < np ^ trunc_log p n n.

Lemma trunc_log_max p k j : 1 < pp ^ j kj trunc_log p k.

pi- parts
Testing for membership in set of prime factors.

Canonical nat_pred_pred := Eval hnf in [predType of nat_pred].

Coercion nat_pred_of_nat (p : nat) : nat_pred := pred1 p.

Section NatPreds.

Variables (n : nat) (pi : nat_pred).

Definition negn : nat_pred := [predC pi].

Definition pnat : pred nat := fun m(m > 0) && all (mem pi) (primes m).

Definition partn := \prod_(0 p < n.+1 | p \in pi) p ^ logn p n.

End NatPreds.

Notation "pi ^'" := (negn pi) (at level 2, format "pi ^'") : nat_scope.

Notation "pi .-nat" := (pnat pi) (at level 2, format "pi .-nat") : nat_scope.

Notation "n `_ pi" := (partn n pi) : nat_scope.

Section PnatTheory.

Implicit Types (n p : nat) (pi rho : nat_pred).

Lemma negnK pi : pi^'^' =i pi.

Lemma eq_negn pi1 pi2 : pi1 =i pi2pi1^' =i pi2^'.

Lemma eq_piP m n : \pi(m) =i \pi(n) \pi(m) = \pi(n).

Lemma part_gt0 pi n : 0 < n`_pi.
Hint Resolve part_gt0.

Lemma sub_in_partn pi1 pi2 n :
  {in \pi(n), {subset pi1 pi2}}n`_pi1 %| n`_pi2.

Lemma eq_in_partn pi1 pi2 n : {in \pi(n), pi1 =i pi2}n`_pi1 = n`_pi2.

Lemma eq_partn pi1 pi2 n : pi1 =i pi2n`_pi1 = n`_pi2.

Lemma partnNK pi n : n`_pi^'^' = n`_pi.

Lemma widen_partn m pi n :
  n mn`_pi = \prod_(0 p < m.+1 | p \in pi) p ^ logn p n.

Lemma partn0 pi : 0`_pi = 1.

Lemma partn1 pi : 1`_pi = 1.

Lemma partnM pi m n : m > 0 → n > 0 → (m × n)`_pi = m`_pi × n`_pi.

Lemma partnX pi m n : (m ^ n)`_pi = m`_pi ^ n.

Lemma partn_dvd pi m n : n > 0 → m %| nm`_pi %| n`_pi.

Lemma p_part p n : n`_p = p ^ logn p n.

Lemma p_part_eq1 p n : (n`_p == 1) = (p \notin \pi(n)).

Lemma p_part_gt1 p n : (n`_p > 1) = (p \in \pi(n)).

Lemma primes_part pi n : primes n`_pi = filter (mem pi) (primes n).

Lemma filter_pi_of n m : n < mfilter \pi(n) (index_iota 0 m) = primes n.

Lemma partn_pi n : n > 0 → n`_\pi(n) = n.

Lemma partnT n : n > 0 → n`_predT = n.

Lemma partnC pi n : n > 0 → n`_pi × n`_pi^' = n.

Lemma dvdn_part pi n : n`_pi %| n.

Lemma logn_part p m : logn p m`_p = logn p m.

Lemma partn_lcm pi m n : m > 0 → n > 0 → (lcmn m n)`_pi = lcmn m`_pi n`_pi.

Lemma partn_gcd pi m n : m > 0 → n > 0 → (gcdn m n)`_pi = gcdn m`_pi n`_pi.

Lemma partn_biglcm (I : finType) (P : pred I) F pi :
    ( i, P iF i > 0) →
  (\big[lcmn/1%N]_(i | P i) F i)`_pi = \big[lcmn/1%N]_(i | P i) (F i)`_pi.

Lemma partn_biggcd (I : finType) (P : pred I) F pi :
    #|SimplPred P| > 0 → ( i, P iF i > 0) →
  (\big[gcdn/0]_(i | P i) F i)`_pi = \big[gcdn/0]_(i | P i) (F i)`_pi.

Lemma sub_in_pnat pi rho n :
  {in \pi(n), {subset pi rho}}pi.-nat nrho.-nat n.

Lemma eq_in_pnat pi rho n : {in \pi(n), pi =i rho}pi.-nat n = rho.-nat n.

Lemma eq_pnat pi rho n : pi =i rhopi.-nat n = rho.-nat n.

Lemma pnatNK pi n : pi^'^'.-nat n = pi.-nat n.

Lemma pnatI pi rho n : [predI pi & rho].-nat n = pi.-nat n && rho.-nat n.

Lemma pnat_mul pi m n : pi.-nat (m × n) = pi.-nat m && pi.-nat n.

Lemma pnat_exp pi m n : pi.-nat (m ^ n) = pi.-nat m || (n == 0).

Lemma part_pnat pi n : pi.-nat n`_pi.

Lemma pnatE pi p : prime ppi.-nat p = (p \in pi).

Lemma pnat_id p : prime pp.-nat p.

Lemma coprime_pi' m n : m > 0 → n > 0 → coprime m n = \pi(m)^'.-nat n.

Lemma pnat_pi n : n > 0 → \pi(n).-nat n.

Lemma pi_of_dvd m n : m %| nn > 0 → {subset \pi(m) \pi(n)}.

Lemma pi_ofM m n : m > 0 → n > 0 → \pi(m × n) =i [predU \pi(m) & \pi(n)].

Lemma pi_of_part pi n : n > 0 → \pi(n`_pi) =i [predI \pi(n) & pi].

Lemma pi_of_exp p n : n > 0 → \pi(p ^ n) = \pi(p).

Lemma pi_of_prime p : prime p\pi(p) =i (p : nat_pred).

Lemma p'natEpi p n : n > 0 → p^'.-nat n = (p \notin \pi(n)).

Lemma p'natE p n : prime pp^'.-nat n = ~~ (p %| n).

Lemma pnatPpi pi n p : pi.-nat np \in \pi(n)p \in pi.

Lemma pnat_dvd m n pi : m %| npi.-nat npi.-nat m.

Lemma pnat_div m n pi : m %| npi.-nat npi.-nat (n %/ m).

Lemma pnat_coprime pi m n : pi.-nat mpi^'.-nat ncoprime m n.

Lemma p'nat_coprime pi m n : pi^'.-nat mpi.-nat ncoprime m n.

Lemma sub_pnat_coprime pi rho m n :
  {subset rho pi^'}pi.-nat mrho.-nat ncoprime m n.

Lemma coprime_partC pi m n : coprime m`_pi n`_pi^'.

Lemma pnat_1 pi n : pi.-nat npi^'.-nat nn = 1.

Lemma part_pnat_id pi n : pi.-nat nn`_pi = n.

Lemma part_p'nat pi n : pi^'.-nat nn`_pi = 1.

Lemma partn_eq1 pi n : n > 0 → (n`_pi == 1) = pi^'.-nat n.

Lemma pnatP pi n :
  n > 0 → reflect ( p, prime pp %| np \in pi) (pi.-nat n).

Lemma pi_pnat pi p n : p.-nat np \in pipi.-nat n.

Lemma p_natP p n : p.-nat n{k | n = p ^ k}.

Lemma pi'_p'nat pi p n : pi^'.-nat np \in pip^'.-nat n.

Lemma pi_p'nat p pi n : pi.-nat np \in pi^'p^'.-nat n.

Lemma partn_part pi rho n : {subset pi rho}n`_rho`_pi = n`_pi.

Lemma partnI pi rho n : n`_[predI pi & rho] = n`_pi`_rho.

Lemma odd_2'nat n : odd n = 2^'.-nat n.

End PnatTheory.
Hint Resolve part_gt0.

Properties of the divisors list.

Lemma divisors_correct n : n > 0 →
  [/\ uniq (divisors n), sorted leq (divisors n)
    & d, (d \in divisors n) = (d %| n)].

Lemma sorted_divisors n : sorted leq (divisors n).

Lemma divisors_uniq n : uniq (divisors n).

Lemma sorted_divisors_ltn n : sorted ltn (divisors n).

Lemma dvdn_divisors d m : 0 < m(d %| m) = (d \in divisors m).

Lemma divisor1 n : 1 \in divisors n.

Lemma divisors_id n : 0 < nn \in divisors n.

Big sum / product lemmas

Lemma dvdn_sum d I r (K : pred I) F :
  ( i, K id %| F i) → d %| \sum_(i <- r | K i) F i.

Lemma dvdn_partP n m : 0 < n
  reflect ( p, p \in \pi(n)n`_p %| m) (n %| m).

Lemma modn_partP n a b : 0 < n
  reflect ( p : nat, p \in \pi(n)a = b %[mod n`_p]) (a == b %[mod n]).

The Euler totient function

Lemma totientE n :
  n > 0 → totient n = \prod_(p <- primes n) (p.-1 × p ^ (logn p n).-1).

Lemma totient_gt0 n : (0 < totient n) = (0 < n).

Lemma totient_pfactor p e :
  prime pe > 0 → totient (p ^ e) = p.-1 × p ^ e.-1.

Lemma totient_coprime m n :
  coprime m ntotient (m × n) = totient m × totient n.

Lemma totient_count_coprime n : totient n = \sum_(0 d < n) coprime n d.