(* ========================================================================= *)
(* Some divisibility properties of certain linear integer recurrences.       *)
(* ========================================================================= *)

needs "Library/prime.ml";;
needs "Library/integer.ml";;

prioritize_int();;

(* ------------------------------------------------------------------------- *)
(* A customized induction principle.                                         *)
(* ------------------------------------------------------------------------- *)

let INDUCT_SPECIAL = 
prove (`!P. (!n. P 0 n) /\ (!m n. P m n <=> P n m) /\ (!m n. P m n ==> P n (m + n)) ==> !m n. P m n`,
GEN_TAC THEN STRIP_TAC THEN REPEAT GEN_TAC THEN WF_INDUCT_TAC `m + n:num` THEN ASM_CASES_TAC `m = 0` THENL [ASM_MESON_TAC[]; ALL_TAC] THEN ASM_CASES_TAC `n = 0` THENL [ASM_MESON_TAC[]; ALL_TAC] THEN DISJ_CASES_THEN MP_TAC (ARITH_RULE `m <= n:num \/ n <= m`) THEN REWRITE_TAC[LE_EXISTS] THEN ONCE_REWRITE_TAC[ADD_SYM] THEN DISCH_THEN(X_CHOOSE_THEN `p:num` SUBST_ALL_TAC) THENL [ALL_TAC; ASM (GEN_REWRITE_TAC I) []] THEN MATCH_MP_TAC(ASSUME `!m n:num. P m n ==> P n (m + n)`) THEN FIRST_X_ASSUM MATCH_MP_TAC THEN ASM_ARITH_TAC);;
(* ------------------------------------------------------------------------- *) (* The main results; to literally apply integer gcd we need nonnegativity. *) (* ------------------------------------------------------------------------- *)
let INT_DIVISORS_RECURRENCE = 
prove (`!G a b. G(0) = &0 /\ G(1) = &1 /\ coprime(a,b) /\ (!n. G(n + 2) = a * G(n + 1) + b * G(n)) ==> !d m n. d divides (G m) /\ d divides (G n) <=> d divides G(gcd(m,n))`,
REPEAT GEN_TAC THEN STRIP_TAC THEN SUBGOAL_THEN `!n. coprime(G(n + 1),b)` ASSUME_TAC THENL [INDUCT_TAC THEN ASM_REWRITE_TAC[ARITH; ARITH_RULE `SUC n + 1 = n + 2`] THEN REPEAT(POP_ASSUM MP_TAC) THEN NUMBER_TAC; ALL_TAC] THEN SUBGOAL_THEN `!n. coprime(G(n + 1),G n)` ASSUME_TAC THENL [INDUCT_TAC THENL [ASM_REWRITE_TAC[ARITH] THEN NUMBER_TAC; ALL_TAC] THEN REPEAT(FIRST_X_ASSUM(ASSUME_TAC o SPEC `n:num`)) THEN ASM_REWRITE_TAC[ADD1; ARITH_RULE `(n + 1) + 1 = n + 2`] THEN REPEAT(POP_ASSUM MP_TAC) THEN INTEGER_TAC; ALL_TAC] THEN SUBGOAL_THEN `!m p. G(m + 1 + p) = G(m + 1) * G(p + 1) + b * G(m) * G(p)` ASSUME_TAC THENL [INDUCT_TAC THENL [ASM_REWRITE_TAC[ADD_CLAUSES; ADD_AC] THEN INTEGER_TAC; ALL_TAC] THEN ASM_REWRITE_TAC[ARITH_RULE `SUC m + 1 + p = (m + p) + 2`] THEN ASM_REWRITE_TAC[ARITH_RULE `SUC m + 1 = m + 2`] THEN ASM_REWRITE_TAC[ARITH_RULE `(m + p) + 1 = m + 1 + p`] THEN INDUCT_TAC THEN ASM_REWRITE_TAC[ARITH; ADD_CLAUSES] THEN ASM_REWRITE_TAC[ARITH_RULE `SUC(m + p) = m + 1 + p`] THEN ASM_REWRITE_TAC[ARITH_RULE `SUC(m + 1) = m + 2`; ARITH] THEN REWRITE_TAC[ADD1] THEN ARITH_TAC; ALL_TAC] THEN SUBGOAL_THEN `!m p:num. gcd(G(m + p),G m) = gcd(G m,G p)` ASSUME_TAC THENL [INDUCT_TAC THEN REWRITE_TAC[ADD_CLAUSES; EQT_INTRO(SPEC_ALL INT_GCD_SYM)] THEN ASM_REWRITE_TAC[ADD1; ARITH_RULE `(m + p) + 1 = m + 1 + p`] THEN GEN_TAC THEN SIMP_TAC[INT_GCD_POS; GSYM INT_DIVIDES_ANTISYM_POS] THEN MP_TAC(SPEC `m:num` (ASSUME `!n. coprime(G(n + 1),b)`)) THEN MP_TAC(SPEC `m:num` (ASSUME `!n. coprime(G(n + 1),G n)`)) THEN INTEGER_TAC; ALL_TAC] THEN GEN_TAC THEN MATCH_MP_TAC INDUCT_SPECIAL THEN REPEAT CONJ_TAC THENL [ASM_REWRITE_TAC[GCD_0; INT_DIVIDES_0]; MESON_TAC[GCD_SYM]; ALL_TAC] THEN ASM_MESON_TAC[GCD_ADD; INT_DIVIDES_GCD; INT_GCD_SYM; ADD_SYM; GCD_SYM]);;
let INT_GCD_RECURRENCE = 
prove (`!G a b. G(0) = &0 /\ G(1) = &1 /\ coprime(a,b) /\ (!n. G(n + 2) = a * G(n + 1) + b * G(n)) /\ (!n. &0 <= G n) ==> !m n. gcd(G m,G n) = G(gcd(m,n))`,
REPEAT GEN_TAC THEN DISCH_TAC THEN ASM_SIMP_TAC[GSYM INT_DIVIDES_ANTISYM_POS; INT_GCD_POS] THEN REWRITE_TAC[INT_DIVIDES_ANTISYM_DIVISORS; INT_DIVIDES_GCD] THEN ASM_MESON_TAC[INT_DIVISORS_RECURRENCE]);;
(* ------------------------------------------------------------------------- *) (* Natural number variants of the same results. *) (* ------------------------------------------------------------------------- *)
let GCD_RECURRENCE = 
prove (`!G a b. G(0) = 0 /\ G(1) = 1 /\ coprime(a,b) /\ (!n. G(n + 2) = a * G(n + 1) + b * G(n)) ==> !m n. gcd(G m,G n) = G(gcd(m,n))`,
REPEAT STRIP_TAC THEN MP_TAC(SPECL [`& o (G:num->num)`; `&a:int`; `&b:int`] INT_GCD_RECURRENCE) THEN ASM_REWRITE_TAC[o_THM; GSYM INT_OF_NUM_ADD; GSYM INT_OF_NUM_MUL] THEN ASM_SIMP_TAC[GSYM num_coprime; INT_POS; GSYM NUM_GCD; INT_OF_NUM_EQ]);;
let DIVISORS_RECURRENCE = 
prove (`!G a b. G(0) = 0 /\ G(1) = 1 /\ coprime(a,b) /\ (!n. G(n + 2) = a * G(n + 1) + b * G(n)) ==> !d m n. d divides (G m) /\ d divides (G n) <=> d divides G(gcd(m,n))`,
REWRITE_TAC[GSYM DIVIDES_GCD] THEN MESON_TAC[DIVISORS_EQ; GCD_RECURRENCE]);;
(* ------------------------------------------------------------------------- *) (* Application 1: Mersenne numbers. *) (* ------------------------------------------------------------------------- *)
let GCD_MERSENNE = 
prove (`!m n. gcd(2 EXP m - 1,2 EXP n - 1) = 2 EXP (gcd(m,n)) - 1`,
SIMP_TAC[GSYM INT_OF_NUM_EQ; NUM_GCD; GSYM INT_OF_NUM_SUB; GSYM INT_OF_NUM_POW; EXP_LT_0; ARITH; ARITH_RULE `1 <= n <=> 0 < n`] THEN MATCH_MP_TAC INT_GCD_RECURRENCE THEN MAP_EVERY EXISTS_TAC [`&3`; `-- &2`] THEN REWRITE_TAC[INT_POW_ADD; INT_LE_SUB_LADD] THEN CONV_TAC INT_REDUCE_CONV THEN REPEAT CONJ_TAC THENL [REWRITE_TAC[GSYM(INT_REDUCE_CONV `&2 * &2 - &1`)] THEN SPEC_TAC(`&2`,`t:int`) THEN INTEGER_TAC; INT_ARITH_TAC; GEN_TAC THEN MATCH_MP_TAC INT_POW_LE_1 THEN INT_ARITH_TAC]);;
let DIVIDES_MERSENNE = 
prove (`!m n. (2 EXP m - 1) divides (2 EXP n - 1) <=> m divides n`,
REPEAT GEN_TAC THEN REWRITE_TAC[DIVIDES_GCD_LEFT; GCD_MERSENNE] THEN SIMP_TAC[EXP_EQ_0; EQ_EXP; ARITH_EQ; ARITH_RULE `~(x = 0) /\ ~(y = 0) ==> (x - 1 = y - 1 <=> x = y)`]);;
(* ------------------------------------------------------------------------- *) (* Application 2: the Fibonacci series. *) (* ------------------------------------------------------------------------- *)
let fib = define
 `fib 0 = 0 /\ fib 1 = 1 /\ !n. fib(n + 2) = fib(n + 1) + fib(n)`;;
let GCD_FIB = 
prove (`!m n. gcd(fib m,fib n) = fib(gcd(m,n))`,
MATCH_MP_TAC GCD_RECURRENCE THEN REPEAT(EXISTS_TAC `1`) THEN REWRITE_TAC[fib; COPRIME_1] THEN ARITH_TAC);;
let FIB_EQ_0 = 
prove (`!n. fib n = 0 <=> n = 0`,
MATCH_MP_TAC num_INDUCTION THEN REWRITE_TAC[fib] THEN MATCH_MP_TAC num_INDUCTION THEN REWRITE_TAC[fib; ARITH_RULE `SUC(SUC n) = n + 2`; ADD_EQ_0] THEN SIMP_TAC[ADD1; ADD_EQ_0; ARITH_EQ] THEN CONV_TAC NUM_REDUCE_CONV THEN REWRITE_TAC[fib; ARITH_EQ]);;
let FIB_INCREASES_LE = 
prove (`!m n. m <= n ==> fib m <= fib n`,
MATCH_MP_TAC TRANSITIVE_STEPWISE_LE THEN REWRITE_TAC[LE_REFL; LE_TRANS] THEN MATCH_MP_TAC num_INDUCTION THEN REWRITE_TAC[fib; ARITH] THEN REWRITE_TAC[ADD1; fib; ARITH_RULE `(n + 1) + 1 = n + 2`] THEN ARITH_TAC);;
let FIB_INCREASES_LT = 
prove (`!m n. 2 <= m /\ m < n ==> fib m < fib n`,
INDUCT_TAC THEN REWRITE_TAC[ARITH] THEN REPEAT STRIP_TAC THEN TRANS_TAC LTE_TRANS `fib(m + 2)` THEN ASM_SIMP_TAC[FIB_INCREASES_LE; ARITH_RULE `m + 2 <= n <=> SUC m < n`] THEN REWRITE_TAC[fib; ADD1; ARITH_RULE `m < m + n <=> ~(n = 0)`; FIB_EQ_0] THEN ASM_ARITH_TAC);;
let FIB_EQ_1 = 
prove (`!n. fib n = 1 <=> n = 1 \/ n = 2`,
MATCH_MP_TAC num_INDUCTION THEN REWRITE_TAC[fib; ARITH] THEN MATCH_MP_TAC num_INDUCTION THEN REWRITE_TAC[fib; ARITH] THEN REWRITE_TAC[fib; ARITH_RULE `SUC(SUC n) = n + 2`] THEN REWRITE_TAC[FIB_EQ_0; ADD_EQ_0; ARITH; ARITH_RULE `m + n = 1 <=> m = 0 /\ n = 1 \/ m = 1 /\ n = 0`] THEN ARITH_TAC);;
let DIVIDES_FIB = 
prove (`!m n. (fib m) divides (fib n) <=> m divides n \/ n = 0 \/ m = 2`,
REPEAT GEN_TAC THEN REWRITE_TAC[DIVIDES_GCD_LEFT; GCD_FIB] THEN MP_TAC(SPECL [`gcd(m:num,n)`; `m:num`] DIVIDES_LE) THEN REWRITE_TAC[GCD] THEN ASM_CASES_TAC `m = 0` THEN ASM_REWRITE_TAC[GCD_0; fib; FIB_EQ_0; ARITH] THEN ASM_CASES_TAC `n = 0` THEN ASM_REWRITE_TAC[GCD_0] THEN ASM_CASES_TAC `gcd(m:num,n) = m` THEN ASM_REWRITE_TAC[LE_LT] THEN ASM_CASES_TAC `gcd(m:num,n) = 0` THENL [ASM_MESON_TAC[GCD_ZERO]; ALL_TAC] THEN ASM_CASES_TAC `m:num = n` THEN ASM_REWRITE_TAC[GCD_REFL; LT_REFL] THEN ASM_CASES_TAC `2 <= gcd(m,n)` THENL [MP_TAC(SPECL [`gcd(m:num,n)`; `m:num`] FIB_INCREASES_LT) THEN ASM_ARITH_TAC; ASM_CASES_TAC `gcd(m,n) = 1` THENL [ASM_REWRITE_TAC[]; ASM_ARITH_TAC] THEN DISCH_TAC THEN CONV_TAC(LAND_CONV SYM_CONV) THEN REWRITE_TAC[FIB_EQ_1; fib] THEN ASM_ARITH_TAC]);;
(* ------------------------------------------------------------------------- *) (* Application 3: solutions of the Pell equation x^2 = (a^2 - 1) y^2 + 1. *) (* All solutions are of the form (pellx a n,pelly a n); see Examples/pell.ml *) (* ------------------------------------------------------------------------- *)
let pellx = define
  `(!a. pellx a 0 = 1) /\
   (!a. pellx a 1 = a) /\
   (!a n. pellx a (n + 2) = 2 * a * pellx a (n + 1) - pellx a n)`;;
let pelly = define
  `(!a. pelly a 0 = 0) /\
   (!a. pelly a 1 = 1) /\
   (!a n. pelly a (n + 2) = 2 * a * pelly a (n + 1) - pelly a (n))`;;
let PELLY_INCREASES = 
prove (`!a n. ~(a = 0) ==> pelly a n <= pelly a (n + 1)`,
GEN_TAC THEN REWRITE_TAC[RIGHT_FORALL_IMP_THM] THEN DISCH_TAC THEN INDUCT_TAC THEN ASM_SIMP_TAC[pelly; ARITH; LE_1; ADD1; ARITH_RULE `(n + 1) + 1 = n + 2`] THEN TRANS_TAC LE_TRANS `2 * pelly a (n + 1) - pelly a n` THEN CONJ_TAC THENL [ASM_ARITH_TAC; ALL_TAC] THEN MATCH_MP_TAC(ARITH_RULE `a:num <= b ==> a - c <= b - c`) THEN REWRITE_TAC[MULT_ASSOC; LE_MULT_RCANCEL] THEN ASM_ARITH_TAC);;
let GCD_PELLY = 
prove (`!a m n. ~(a = 0) ==> gcd(pelly a m,pelly a n) = pelly a (gcd(m,n))`,
GEN_TAC THEN REWRITE_TAC[RIGHT_FORALL_IMP_THM] THEN DISCH_TAC THEN REWRITE_TAC[GSYM INT_OF_NUM_EQ; NUM_GCD] THEN MATCH_MP_TAC INT_GCD_RECURRENCE THEN MAP_EVERY EXISTS_TAC [`&2 * &a:int`; `-- &1:int`] THEN REWRITE_TAC[pelly; INT_POS; INT_COPRIME_NEG; INT_COPRIME_1] THEN GEN_TAC THEN REWRITE_TAC[INT_OF_NUM_MUL; MULT_ASSOC] THEN REWRITE_TAC[INT_ARITH `a + -- &1 * b:int = a - b`] THEN MATCH_MP_TAC(GSYM INT_OF_NUM_SUB) THEN TRANS_TAC LE_TRANS `1 * pelly a (n + 1)` THEN REWRITE_TAC[LE_MULT_RCANCEL] THEN ASM_SIMP_TAC[MULT_CLAUSES; PELLY_INCREASES] THEN ASM_ARITH_TAC);;