Library prime

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.      
         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 [::].         
          phi 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.      


 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 : nat) {struct r} :=
  if r is r'.+2 then edivn2 q.+1 r' else (q, r).

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

Fixpoint elogn2 (e q r : nat) {struct q} :=
  match q, r with
  | 0, _ | _, 0 => (e, q)
  | q'.+1, 1 => elogn2 e.+1 q' q'
  | q'.+1, r'.+2 => elogn2 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 : forall 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 : forall 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 : nat) {struct a} :=
  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_phi_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 phi n := foldr add_phi_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
  forall n, n > 0 -> pd_ok 1 n (prime_decomp n).

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

Lemma primeP : forall p,
  reflect (p > 1 /\ forall d, d %| p -> xpred2 1 p d) (prime p).

Lemma prime_nt_dvdP : forall d p, prime p -> d != 1 -> reflect (d = p) (d %| p).

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

Lemma prime_gt1 : forall p, prime p -> 1 < p.

Lemma prime_gt0 : forall p, prime p -> 0 < p.

Hint Resolve prime_gt1 prime_gt0.

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

Lemma even_prime : forall p, prime p -> p = 2 \/ odd p.

Lemma prime_oddPn : forall p, prime p -> reflect (p = 2) (~~ odd p).

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

Lemma prime_coprime : forall p m, prime p -> coprime p m = ~~ (p %| m).

Lemma dvdn_prime2 : forall p q, prime p -> prime q -> (p %| q) = (p == q).

Lemma euclid : forall m n p, prime p -> (p %| m * n) = (p %| m) || (p %| n).

Lemma euclid1 : forall p, prime p -> (p %| 1) = false.

Lemma euclid_exp : forall m n p, prime p -> (p %| m ^ n) = (p %| m) && (n > 0).

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

Lemma sorted_primes : forall n, sorted ltn (primes n).

Lemma eq_primes : forall m n, (primes m =i primes n) <-> (primes m = primes n).

Lemma primes_uniq : forall n, uniq (primes n).

 The smallest prime divisor 

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

Lemma pdiv_prime : forall n, 1 < n -> prime (pdiv n).

Lemma pdiv_dvd : forall n, pdiv n %| n.

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

Lemma max_pdiv_prime : forall n, n > 1 -> prime (max_pdiv n).

Lemma max_pdiv_dvd : forall n, max_pdiv n %| n.

Lemma pdiv_leq : forall n, 0 < n -> pdiv n <= n.

Lemma max_pdiv_leq : forall n, 0 < n -> max_pdiv n <= n.

Lemma pdiv_gt0 : forall n, 0 < pdiv n.

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

Lemma pdiv_min_dvd : forall m d, 1 < d -> d %| m -> pdiv m <= d.

Lemma max_pdiv_max : forall n p, p \in \pi(n) -> p <= max_pdiv n.

Lemma ltn_pdiv2_prime : forall n, 0 < n -> n < pdiv n ^ 2 -> prime n.

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

Implicit Arguments primePns [n].

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

Lemma primes_mul : forall m n p, m > 0 -> n > 0 ->
  (p \in primes (m * n)) = (p \in primes m) || (p \in primes n).

Lemma primes_exp : forall m n, n > 0 -> primes (m ^ n) = primes m.

Lemma primes_prime : forall p, prime p -> primes p = [::p].

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

Lemma pdiv_id : forall p, prime p -> pdiv p = p.

Lemma pdiv_pfactor : forall p k, prime p -> pdiv (p ^ k.+1) = p.

 "prime" logarithms and p-parts. 

Fixpoint logn_rec (d m r : nat) {struct r} : nat :=
  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 : forall p m,
  logn p m = if [&& prime p, 0 < m & p %| m] then (logn p (m %/ p)).+1 else 0.

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

Lemma ltn_log0 : forall p n, n < p -> logn p n = 0.

Lemma logn0 : forall p, logn p 0 = 0.

Lemma logn1 : forall p, logn p 1 = 0.

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

Lemma pfactor_dvdn : forall p n m,
  prime p -> m > 0 -> (p ^ n %| m) = (n <= logn p m).

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

Lemma logn_prime : forall p q, prime q -> logn p q = (p == q).

Lemma pfactor_coprime : forall p n, prime p -> n > 0 ->
  {m | coprime p m & n = m * p ^ logn p n}.

Lemma pfactorK : forall p n, prime p -> logn p (p ^ n) = n.

Lemma pfactorKpdiv : forall p n,
  prime p -> logn (pdiv (p ^ n)) (p ^ n) = n.

Lemma dvdn_leq_log : forall p m n, 0 < n -> m %| n -> logn p m <= logn p n.

Lemma logn_gauss : forall p m n, coprime p m -> logn p (m * n) = logn p n.

Lemma logn_mul : forall p m n,
  0 < m -> 0 < n -> logn p (m * n) = logn p m + logn p n.

Lemma logn_exp : forall p m n, logn p (m ^ n) = n * logn p m.

Lemma logn_div : forall p m n, m %| n -> logn p (n %/ m) = logn p n - logn p m.

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

Lemma prime_decompE : forall n,
  prime_decomp n = map (fun p => (p, logn p n)) (primes n).

 pi- parts 

 Testing for membership in set of prime factors. 

Canonical Structure 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.

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

Lemma eq_negn : forall pi1 pi2, pi1 =i pi2 -> pi1^' =i pi2^'.

Lemma eq_piP : forall m n, \pi(m) =i \pi(n) <-> \pi(m) = \pi(n).

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

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

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

Lemma eq_partn : forall (pi1 pi2 : nat_pred) n, pi1 =i pi2 -> n`_pi1 = n`_pi2.

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

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

Lemma partn0 : forall pi, 0`_pi = 1.

Lemma partn1 : forall pi, 1`_pi = 1.

Lemma partn_mul : forall pi m n, m > 0 -> n > 0 -> (m * n)`_pi = m`_pi * n`_pi.

Lemma partn_exp : forall pi m n, (m ^ n)`_pi = m`_pi ^ n.

Lemma partn_dvd : forall pi m n, n > 0 -> m %| n -> m`_pi %| n`_pi.

Lemma p_part : forall p n : nat, n`_p = p ^ logn p n.

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

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

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

Lemma filter_pi_of : forall n m,
  n < m -> filter \pi(n) (index_iota 0 m) = primes n.

Lemma partn_pi : forall n, n > 0 -> n`_\pi(n) = n.

Lemma partnT : forall n, n > 0 -> n`_predT = n.

Lemma partnC : forall pi n, n > 0 -> n`_pi * n`_pi^' = n.

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

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

Lemma partn_lcm : forall pi m n,
  m > 0 -> n > 0 -> (lcmn m n)`_pi = lcmn m`_pi n`_pi.

Lemma partn_gcd : forall pi m n,
  m > 0 -> n > 0 -> (gcdn m n)`_pi = gcdn m`_pi n`_pi.

Lemma partn_biglcm : forall (I : finType) (P : pred I) F pi,
    (forall 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 : forall (I : finType) (P : pred I) F pi,
    #|SimplPred P| > 0 -> (forall i, P i -> F i > 0) ->
  (\big[gcdn/0]_(i | P i) F i)`_pi = \big[gcdn/0]_(i | P i) (F i)`_pi.

Section PiNat.

Implicit Types p n m : nat.
Implicit Type pi rho : nat_pred.

Lemma sub_in_pnat : forall pi rho n,
  {in \pi(n), {subset pi <= rho}} -> pi.-nat n -> rho.-nat n.

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

Lemma eq_pnat : forall pi rho n, pi =i rho -> pi.-nat n = rho.-nat n.

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

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

Lemma pnat_mul : forall pi m n, pi.-nat (m * n) = pi.-nat m && pi.-nat n.

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

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

Lemma pnatE : forall pi p, prime p -> pi.-nat p = (p \in pi).

Lemma pnat_id : forall p, prime p -> p.-nat p.

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

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

Lemma pi_of_dvd : forall m n, m %| n -> n > 0 -> {subset \pi(m) <= \pi(n)}.

Lemma pi_of_muln : forall m n,
  m > 0 -> n > 0 -> \pi(m * n) =i [predU \pi(m) & \pi(n)].

Lemma pi_of_partn : forall pi n, n > 0 -> \pi(n`_pi) =i [predI \pi(n) & pi].

Lemma pi_of_exp : forall p n, n > 0 -> \pi(p ^ n) = \pi(p).

Lemma pi_of_prime : forall p, prime p -> \pi(p) =i (p : nat_pred).

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

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

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

Lemma pnat_dvd : forall m n pi, m %| n -> pi.-nat n -> pi.-nat m.

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

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

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

Lemma sub_pnat_coprime : forall pi rho m n,
  {subset rho <= pi^'} -> pi.-nat m -> rho.-nat n -> coprime m n.

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

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

Lemma part_pnat_id : forall pi n, pi.-nat n -> n`_pi = n.

Lemma part_p'nat : forall pi n, pi^'.-nat n -> n`_pi = 1.

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

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

Lemma pi_pnat : forall pi p n, p.-nat n -> p \in pi -> pi.-nat n.

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

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

Lemma pi_p'nat : forall pi p n, pi.-nat n -> p \in pi^' -> p^'.-nat n.

Lemma partn_part : forall pi rho n,
  {subset pi <= rho} -> n`_rho`_pi = n`_pi.

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

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

End PiNat.

 Properties of the divisors list. 

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

Lemma sorted_divisors : forall n, sorted leq (divisors n).

Lemma divisors_uniq : forall n, uniq (divisors n).

Lemma sorted_divisors_ltn : forall n, sorted ltn (divisors n).

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

Lemma divisor1 : forall n, 1 \in divisors n.

Lemma divisorn : forall n, 0 < n -> n \in divisors n.

 Big sum / product lemmas

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

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

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

 The Euler phi function 

Lemma phiE : forall n,
  n > 0 -> phi n = \prod_(p <- primes n) (p.-1 * p ^ (logn p n).-1).

Lemma phi_gt0 : forall n, (0 < phi n) = (0 < n).

Lemma phi_pfactor : forall p e,
  prime p -> e > 0 -> phi (p ^ e) = p.-1 * p ^ e.-1.

Lemma phi_coprime : forall m n,
  coprime m n -> phi (m * n) = phi m * phi n.

Lemma phi_count_coprime : forall n, phi n = \sum_(0 <= d < n | coprime n d) 1.