# Library mathcomp.ssreflect.prime

(* (c) Copyright 2006-2015 Microsoft Corporation and Inria.
Distributed under the terms of CeCILL-B.                                  *)

Require Import mathcomp.ssreflect.ssreflect.

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. npi == 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 np.
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).

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 × nat Type :=
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) : T Type :=
| 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.

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

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 %| p xpred2 1 p d) (prime p).

Lemma prime_nt_dvdP d p : prime p d != 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 p p = 2 odd p.

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

Lemma odd_prime_gt2 p : odd p prime p p > 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 p coprime p m = ~~ (p %| m).

Lemma dvdn_prime2 p q : prime p prime 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 < n prime (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 < n pdiv n n.

Lemma max_pdiv_leq n : 0 < n max_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 < d d %| m pdiv m d.

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

Lemma ltn_pdiv2_prime n : 0 < n n < 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 p primes 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 p pdiv p = p.

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

Primes are unbounded.

Lemma dvdn_fact m n : 0 < m n m %| n!.

Lemma prime_above m : {p | m < p & prime 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 < p logn 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 p m > 0 (p ^ n %| m) = (n logn p m).

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

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

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

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

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

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

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

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

Lemma lognM p m n : 0 < m 0 < n logn 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 %| n logn 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 p logn 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 < n let k := trunc_log p n in p ^ k n < p ^ k.+1.

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

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

Lemma trunc_log_max p k j : 1 < p p ^ j k j 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 pi2 pi1^' =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 pi2 n_pi1 = n_pi2.

Lemma partnNK pi n : n_pi^'^' = n_pi.

Lemma widen_partn m pi n :
n m n_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 %| n m_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 < m filter \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 i F 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 i F 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 n rho.-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 rho pi.-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 p pi.-nat p = (p \in pi).

Lemma pnat_id p : prime p p.-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 %| n n > 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 p p^'.-nat n = ~~ (p %| n).

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

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

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

Lemma pnat_coprime pi m n : pi.-nat m pi^'.-nat n coprime m n.

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

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

Lemma coprime_partC pi m n : coprime m_pi n_pi^'.

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

Lemma part_pnat_id pi n : pi.-nat n n_pi = n.

Lemma part_p'nat pi n : pi^'.-nat n n_pi = 1.

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

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

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

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

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

Lemma pi_p'nat p pi n : pi.-nat n p \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 < n n \in divisors n.

Big sum / product lemmas

Lemma dvdn_sum d I r (K : pred I) F :
( i, K i d %| 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