# Library mathcomp.ssreflect.prime

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

From mathcomp Require Import ssreflect ssrfun ssrbool eqtype ssrnat seq path.
From mathcomp Require Import fintype 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 p e == the value p ^ e of a prime factor (p, e). NumFactor f == print version of a prime factor, converting the prime component to a Num (which can print large values). 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 <= 1 or m is 0. up_log p m == the smallest e such that m <= p ^ e, or 0 if p <= 1 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. The type of A must coerce to finpred_sort (e.g., by coercing to {set T}) and not merely implement the predType interface (as seq T does).
• > The expression #|A| will only appear in \pi(A) after simplification collapses the coercion, 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 arithmetic 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).

Module Import PrimeDecompAux.

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.
Arguments elogn2 : simpl nomatch.

Variant 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.

Variant 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).

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

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

Import NatTrec.

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

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

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

End PrimeDecompAux.

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

Definition pfactor p e := p ^ e.

Section prime_decomp.

Import NatTrec.

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].

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_arg := nat.
Coercion pi_arg_of_nat (n : nat) : pi_arg := n.
Coercion pi_arg_of_fin_pred T pT (A : @fin_pred_sort T pT) : pi_arg := #|A|.
Arguments pi_arg_of_nat n /.
Arguments pi_arg_of_fin_pred {T pT} A /.
Definition pi_of (n : pi_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).

Arguments primeP {p}.
Arguments primePn {n}.

Lemma prime_gt1 p : prime p 1 < p.

Lemma prime_gt0 p : prime p 0 < p.

#[global] Hint Resolve prime_gt1 prime_gt0 : core.

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_dvd1 p : prime p (p %| 1) = false.

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

Lemma Euclid_dvd_prod (I : Type) (r : seq I) (P : pred I) (f : I nat) p :
prime p
p %| \prod_(i <- r | P i) f i = \big[orb/false]_(i <- r | P i) (p %| f i).

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.
#[global] Hint Resolve pdiv_gt0 max_pdiv_gt0 : core.

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).

Arguments primePns {n}.

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

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

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

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

Lemma coprime_has_primes m n :
0 < m 0 < n coprime m n = ~~ has [in 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 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.
#[global] Hint Resolve pfactor_gt0 : core.

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 logn_coprime p m : coprime p m logn p m = 0.

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.
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 if p 1 then 0 else loop n n.

Lemma trunc_log0 p : trunc_log p 0 = 0.

Lemma trunc_log1 p : trunc_log p 1 = 0.

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_logP p n : 1 < p 0 < n p ^ trunc_log p n n.

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

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

Lemma trunc_log_eq0 p n : (trunc_log p n == 0) = (p 1) || (n p.-1).

Lemma trunc_log_gt0 p n : (0 < trunc_log p n) = (1 < p) && (p.-1 < n).

Lemma trunc_log0n n : trunc_log 0 n = 0.

Lemma trunc_log1n n : trunc_log 1 n = 0.

Lemma leq_trunc_log p m n : m n trunc_log p m trunc_log p n.

Lemma trunc_log_eq p n k : 1 < p p ^ n k < p ^ n.+1 trunc_log p k = n.

Lemma trunc_lognn p : 1 < p trunc_log p p = 1.

Lemma trunc_expnK p n : 1 < p trunc_log p (p ^ n) = n.

Lemma trunc_logMp p n : 1 < p 0 < n
trunc_log p (p × n) = (trunc_log p n).+1.

Lemma trunc_log2_double n : 0 < n trunc_log 2 n.*2 = (trunc_log 2 n).+1.

Lemma trunc_log2S n : 1 < n trunc_log 2 n = (trunc_log 2 n./2).+1.

Truncated up real logarithm

Definition up_log p n :=
if (p 1) then 0 else
let v := trunc_log p n in if n p ^ v then v else v.+1.

Lemma up_log0 p : up_log p 0 = 0.

Lemma up_log1 p : up_log p 1 = 0.

Lemma up_log_eq0 p n : (up_log p n == 0) = (p 1) || (n 1).

Lemma up_log_gt0 p n : (0 < up_log p n) = (1 < p) && (1 < n).

Lemma up_log_bounds p n :
1 < p 1 < n let k := up_log p n in p ^ k.-1 < n p ^ k.

Lemma up_logP p n : 1 < p n p ^ up_log p n.

Lemma up_log_gtn p n : 1 < p 1 < n p ^ (up_log p n).-1 < n.

Lemma up_log_min p k j : 1 < p k p ^ j up_log p k j.

Lemma leq_up_log p m n : m n up_log p m up_log p n.

Lemma up_log_eq p n k : 1 < p p ^ n < k p ^ n.+1 up_log p k = n.+1.

Lemma up_lognn p : 1 < p up_log p p = 1.

Lemma up_expnK p n : 1 < p up_log p (p ^ n) = n.

Lemma up_logMp p n : 1 < p 0 < n up_log p (p × n) = (up_log p n).+1.

Lemma up_log2_double n : 0 < n up_log 2 n.*2 = (up_log 2 n).+1.

Lemma up_log2S n : 0 < n up_log 2 n.+1 = (up_log 2 (n./2.+1)).+1.

Lemma up_log_trunc_log p n :
1 < p 1 < n up_log p n = (trunc_log p n.-1).+1.

Lemma trunc_log_up_log p n :
1 < p 0 < n trunc_log p n = (up_log p n.+1).-1.

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 [in 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 : core.

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 eq_partn_from_log m n (pi : nat_pred) : 0 < m 0 < n
{in pi, logn^~ m =1 logn^~ n} m_pi = n_pi.

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 [in 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 eqn_from_log m n : 0 < m 0 < n logn^~ m =1 logn^~ n m = 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 logn_gcd p m n : 0 < m 0 < n
logn p (gcdn m n) = minn (logn p m) (logn p n).

Lemma logn_lcm p m n : 0 < m 0 < n
logn p (lcmn m n) = maxn (logn p m) (logn p n).

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 pnatM pi m n : pi.-nat (m × n) = pi.-nat m && pi.-nat n.

Lemma pnatX 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.
#[global] Hint Resolve part_gt0 : core.

Properties of the divisors list.
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