@@ -2,12 +2,25 @@ From Stdlib Require Import String NArith.
22From mathcomp Require Import all_ssreflect.
33Require Import digitn.
44
5-
65(***************************************************************************** *)
76(* *)
7+ (* sum_fact b n == the sum of the factorial of its digits in base b *)
8+ (* *)
89(* factorion n == n is equal to the sum of the factorial of its digits in *)
9- (* base 10 *)
10- (* This file shows that there is only 4 of these 1, 2, 145, 40585 *)
10+ (* base 10 (i.e n = sum_fact 10 n) *)
11+ (* *)
12+ (* Factorions are fixed-points of the function sum_fact 10 (1-cycle). *)
13+ (* This file shows that there are only 4 of these 1, 2, 145, 40585. *)
14+ (* *)
15+ (* Amical factorions correspond to 2-cycles, we show that 781 and 782 are *)
16+ (* amical factorions. *)
17+ (* *)
18+ (* Magical factorions correspond to 3-cycles, we show that 169 is a magical *)
19+ (* factorions. *)
20+ (* *)
21+ (* Finally we show that these are the only cycles since when iterating the *)
22+ (* application of sum_fact 10 from any n, we show that we always encounter *)
23+ (* one of these numbers [:: 1, 2, 145, 40585, 781, 782, 169] *)
1124(* *)
1225(***************************************************************************** *)
1326
@@ -21,6 +34,17 @@ Let v40585 := Eval vm_compute in 5 + 10 * 4058.
2134
2235Definition sum_fact b (n : nat) := \sum_(i < ndigits b n) (digitn b n i) `!.
2336
37+ Lemma sum_fact0 b : sum_fact b 0 = 1.
38+ Proof .
39+ by rewrite /sum_fact ndigits0 big_ord_recl big_ord0 // digitn0 mod0n.
40+ Qed .
41+
42+ Lemma sum_fact_small b d : 1 < b -> d < b -> sum_fact b d = d `!.
43+ Proof .
44+ move=> b_gt1 dLb; rewrite /sum_fact ndigits_small // big_ord_recl big_ord0.
45+ by rewrite digitn0 modn_small // addn0.
46+ Qed .
47+
2448Lemma sum_factMD b n d :
2549 0 < n -> 1 < b -> d < b -> sum_fact b (b * n + d) = d `! + sum_fact b n.
2650Proof .
@@ -383,6 +407,234 @@ case.
383407by [].
384408Time Qed .
385409
410+ Definition Nfact_small d :=
411+ match d with
412+ | 9%num => 362880%num
413+ | 8%num => 40320%num
414+ | 7%num => 5040%num
415+ | 6%num => 720%num
416+ | 5%num => 120%num
417+ | 4%num => 24%num
418+ | 3%num => 6%num
419+ | 2%num => 2%num
420+ | _ => 1%num
421+ end .
422+
423+ Lemma Nfact_small_spec d :
424+ d < 10 -> Nfact_small (N.of_nat d) = N.of_nat (d `!).
425+ Proof .
426+ case: d => [//|] [//|] [//|] [//|] [//|] [//|] [//|] [//|] [//|] [_|//].
427+ by rewrite 3!factS 3!Nat2N.inj_mul.
428+ Qed .
429+
430+ Fixpoint Nsum_fact10_aux (n : positive) p :=
431+ let p1 := (p / 10)%num in
432+ let r := (p mod 10)%num in
433+ (Nfact_small r + if (p1 =? 0)%num then 0 else
434+ match n with xH => 0%num | xI n1 | xO n1 => Nsum_fact10_aux n1 p1 end )%num.
435+
436+ Lemma NatDiv_div a b : PeanoNat.Nat.div a b = a %/ b.
437+ Proof .
438+ case: b => [|b]; first by rewrite divn0 PeanoNat.Nat.div_0_r.
439+ have [k aLk]:= ubnP a; elim: k a aLk => [[]//|k IH a aLk].
440+ have [bLa|aLb] := leqP b.+1 a; last first.
441+ rewrite PeanoNat.Nat.div_small; last by apply/ltP.
442+ by rewrite divn_small.
443+ rewrite -(subnK bLa) -{2 5}[b.+1]mul1n addnC PeanoNat.Nat.div_add_l // IH //.
444+ by rewrite divnMDl.
445+ rewrite ltnS in aLk.
446+ by rewrite (leq_trans _ aLk) // ltn_subLR // addSnnS leq_addl.
447+ Qed .
448+
449+ Lemma Natmod_mod a b : PeanoNat.Nat.modulo a b = a %% b.
450+ Proof .
451+ case: b => [|b]; first by rewrite modn0 PeanoNat.Nat.mod_0_r.
452+ have [k aLk]:= ubnP a; elim: k a aLk => [[]//|k IH a aLk].
453+ have [bLa|aLb] := leqP b.+1 a; last first.
454+ rewrite PeanoNat.Nat.mod_small; last by apply/ltP.
455+ by rewrite modn_small.
456+ rewrite -(subnK bLa) -{2 5}[b.+1]mul1n PeanoNat.Nat.Div0.mod_add IH //.
457+ by rewrite addnC modnMDl.
458+ rewrite ltnS in aLk.
459+ by rewrite (leq_trans _ aLk) // ltn_subLR // addSnnS leq_addl.
460+ Qed .
461+
462+ Lemma Nsum_fact10_aux_spec n p :
463+ 0 < N.to_nat p <= N.to_nat (Npos n) ->
464+ Nsum_fact10_aux n p = N.of_nat (sum_fact 10 (N.to_nat p)).
465+ Proof .
466+ elim: n p => [n IH p| n IH p| p] pB /=.
467+ - have -> : N.to_nat p = 10 * (N.to_nat (p / 10)) + N.to_nat (p mod 10).
468+ rewrite [LHS](divn_eq _ 10) // mulnC N2Nat.inj_div NatDiv_div.
469+ by rewrite N2Nat.inj_mod Natmod_mod.
470+ have pm10L10 : N.to_nat (p mod 10) < 10.
471+ by rewrite N2Nat.inj_mod Natmod_mod ltn_mod.
472+ case: N.eqb_spec => [->|HNE].
473+ rewrite muln0 add0n N.add_0_r -[(_ mod 10)%num in LHS]N2Nat.id.
474+ by rewrite Nfact_small_spec // sum_fact_small.
475+ have pD10_pos : 0 < N.to_nat (p / 10).
476+ by rewrite -[X in X <> _]N2Nat.id in HNE; case: N.to_nat HNE.
477+ rewrite sum_factMD // Nat2N.inj_add -IH //.
478+ by rewrite -Nfact_small_spec // N2Nat.id.
479+ rewrite pD10_pos N2Nat.inj_div NatDiv_div.
480+ rewrite -[X in _ <= X](mulnK _ (isT : 0 < 10)).
481+ apply: leq_div2r.
482+ apply: leq_trans (_ : N.to_nat (N.pos n~1) <= _); first by case/andP: pB.
483+ suff -> : N.to_nat (N.pos n~1) = (N.to_nat (N.pos n)).*2.+1.
484+ rewrite -muln2 ltn_mul2l // andbT.
485+ by apply/leP/(Pnat.Pos2Nat.inj_le 1)/Pos.le_1_l.
486+ by rewrite /= Pnat.Pos2Nat.inj_xI -mul2n.
487+ - have -> : N.to_nat p = 10 * (N.to_nat (p / 10)) + N.to_nat (p mod 10).
488+ rewrite [LHS](divn_eq _ 10) // mulnC N2Nat.inj_div NatDiv_div.
489+ by rewrite N2Nat.inj_mod Natmod_mod.
490+ have pm10L10 : N.to_nat (p mod 10) < 10.
491+ by rewrite N2Nat.inj_mod Natmod_mod ltn_mod.
492+ case: N.eqb_spec => [->|HNE].
493+ rewrite muln0 add0n N.add_0_r -[(_ mod 10)%num in LHS]N2Nat.id.
494+ by rewrite Nfact_small_spec // sum_fact_small.
495+ have pD10_pos : 0 < N.to_nat (p / 10).
496+ by rewrite -[X in X <> _]N2Nat.id in HNE; case: N.to_nat HNE.
497+ rewrite sum_factMD // Nat2N.inj_add -IH //.
498+ by rewrite -Nfact_small_spec // N2Nat.id.
499+ rewrite pD10_pos N2Nat.inj_div NatDiv_div.
500+ rewrite -[X in _ <= X](mulnK _ (isT : 0 < 10)).
501+ apply: leq_div2r.
502+ apply: leq_trans (_ : N.to_nat (N.pos n~0) <= _); first by case/andP: pB.
503+ suff -> : N.to_nat (N.pos n~0) = (N.to_nat (N.pos n)).*2.
504+ by rewrite -muln2 leq_mul2l // orbT.
505+ by rewrite /= Pnat.Pos2Nat.inj_xO -mul2n.
506+ suff -> : p = 1%num by rewrite /= sum_fact_small.
507+ by apply: N2Nat.inj; case: (N.to_nat p) pB => // [] [].
508+ Qed .
509+
510+ Definition Nsum_fact10 p := if p is Npos n then Nsum_fact10_aux n p else 1%num.
511+
512+ Lemma Nsum_fact10_spec n :
513+ Nsum_fact10 n = N.of_nat (sum_fact 10 (N.to_nat n)).
514+ Proof .
515+ case: n => /= [|p]; first by rewrite sum_fact0.
516+ apply: Nsum_fact10_aux_spec.
517+ rewrite leqnn andbT.
518+ by apply/leP/(Pnat.Pos2Nat.inj_le 1)/Pos.le_1_l.
519+ Qed .
520+
521+ (* The amical factoriom *)
522+
523+ Lemma afactorion871 : sum_fact 10 (sum_fact 10 871) = 871.
524+ Proof .
525+ apply: Nat2N.inj.
526+ rewrite -[ (sum_fact 10 871)]Nat2N.id -Nsum_fact10_spec.
527+ by rewrite -[871 in LHS]Nat2N.id -Nsum_fact10_spec.
528+ Qed .
529+
530+ Lemma afactorion872 : sum_fact 10 (sum_fact 10 872) = 872.
531+ Proof .
532+ apply: Nat2N.inj.
533+ rewrite -[ (sum_fact 10 872)]Nat2N.id -Nsum_fact10_spec.
534+ by rewrite -[872 in LHS]Nat2N.id -Nsum_fact10_spec.
535+ Qed .
536+
537+ (* The magic factorion *)
538+
539+ Lemma mfactorion169 : sum_fact 10 (sum_fact 10 (sum_fact 10 169)) = 169.
540+ Proof .
541+ apply: Nat2N.inj.
542+ rewrite -[sum_fact _ (sum_fact 10 169)]Nat2N.id -Nsum_fact10_spec.
543+ rewrite -[sum_fact 10 169]Nat2N.id -Nsum_fact10_spec.
544+ by rewrite -[169 in LHS]Nat2N.id -Nsum_fact10_spec.
545+ Qed .
546+
547+ (* *)
548+ Lemma sum_fact10_iter_increase n :
549+ exists k, let v := iter k (sum_fact 10) n in v <= (sum_fact 10 v).
550+ Proof .
551+ have [k nLk]:= ubnP n; elim: k n nLk => [[]//|k IH n nLk].
552+ have [nLs|sLn] := leqP n (sum_fact 10 n); first by exists 0.
553+ have /IH[k1 Hk1] // : sum_fact 10 n < k by rewrite -ltnS (leq_trans _ nLk).
554+ by exists k1.+1; rewrite iterSr.
555+ Qed .
556+
557+ Definition in_loop n :=
558+ match n with
559+ 1%num | 2%num | 145%num | 40585%num | 169%num | (* 363601%num | 1454%num | *)
560+ 871%num | (* 45361%num | *) 872%num (* | 45362%num *) => true | _ => false end .
561+
562+ Lemma in_loop_spec n :
563+ in_loop n = (n \in [:: 1; 2; 145; 40585; 871; 872; 169]%num).
564+ Proof . by case: n => //; do 16 (case => //=). Qed .
565+
566+ Fixpoint check_loop p n :=
567+ if in_loop p then true else
568+ if n is n1.+1 then check_loop (Nsum_fact10 p) n1 else false.
569+
570+ Lemma check_loop_spec p n :
571+ check_loop p n -> exists k,
572+ iter k (sum_fact 10) (N.to_nat p) \in [:: 1; 2; 145; v40585; 871; 872; 169].
573+ Proof .
574+ elim: n p => /= [p| n IH p].
575+ rewrite in_loop_spec.
576+ by case: (boolP (_ \in _)) => // /or4P[|||/or4P[|||/orP[]]] // /eqP-> _;
577+ exists 0; rewrite !inE.
578+ rewrite in_loop_spec.
579+ have [|_ /IH [k Hk]] := boolP (_ \in _).
580+ move=> Hp; exists 0; rewrite !inE /=.
581+ by have /or4P[|||/or4P[|||/orP[]]] // := Hp => /eqP->.
582+ rewrite Nsum_fact10_spec Nat2N.id in Hk.
583+ by exists k.+1; rewrite iterSr.
584+ Qed .
585+
586+ Definition cend_in_cycle n m :=
587+ if (n <? m)%num then ~~ check_loop n 60 else false.
588+
589+
590+ Lemma cend_in_cycleN_spec n m :
591+ ~~ cend_in_cycle n m -> (n <? m)%num ->
592+ exists k,
593+ iter k (sum_fact 10) (N.to_nat n)\in [:: 1; 2; 145; v40585; 871; 872; 169].
594+ Proof .
595+ by rewrite /cend_in_cycle; case: (_ <? _)%num;
596+ rewrite ?negbK // => /check_loop_spec.
597+ Qed .
598+
599+ Lemma leq_of_nat n m : (N.of_nat n <= N.of_nat m)%num -> (n <= m).
600+ Proof .
601+ rewrite -{2}(Nat2N.id n) -{2}(Nat2N.id m).
602+ case: N.of_nat => //= p; case: N.of_nat => //= q pLq.
603+ by apply/leP/(Pnat.Pos2Nat.inj_le).
604+ Qed .
605+
606+ (* There is nothing more, we always reach either factorion, amical, or magic *)
607+ Lemma factorion_no_large_cycle n :
608+ exists k, iter k (sum_fact 10) n \in [:: 1; 2; 145; v40585; 871; 872; 169].
609+ Proof .
610+ case: (sum_fact10_iter_increase n) => k.
611+ set v := iter k (sum_fact 10) n => /=.
612+ rewrite leq_eqVlt => /orP[/eqP vE| vLs].
613+ exists k; rewrite -/v.
614+ suff : v \in [:: 1; 2; 145; v40585].
615+ by rewrite !inE => /or4P[] /eqP->; rewrite eqxx.
616+ by rewrite -factorionE; apply/eqP.
617+ suff [k1 Hk1] : exists k,
618+ iter k (sum_fact 10) v \in [:: 1; 2; 145; v40585; 871; 872; 169].
619+ by exists (k1 + k); rewrite iterD.
620+ have [|v_pos] := leqP v 0.
621+ by rewrite leqn0 => /eqP->; exists 1; rewrite /= sum_fact0 !inE.
622+ have nL7 : ndigits 10 v <= 7.
623+ by case: leqP => // /sum_fact10_increasing; rewrite ltnNge leq_eqVlt vLs orbT.
624+ pose get_list := fact_look_up cend_in_cycle.
625+ suff g_nil : get_list (ndigits 10 v).-1 = [::].
626+ have Hv := @fact_look_up_spec cend_in_cycle _ v_pos.
627+ rewrite -/get_list g_nil in_nil in Hv.
628+ have := @cend_in_cycleN_spec (N.of_nat v) (N.of_nat (sum_fact 10 v)).
629+ rewrite Nat2N.id; apply => //; first by case: cend_in_cycle Hv => //; apply.
630+ case: N.ltb_spec => // Hf.
631+ move: vLs; rewrite ltnNge => /negP[].
632+ by apply: leq_of_nat.
633+ case: ndigits nL7 => [_|].
634+ by vm_cast_no_check (refl_equal ([::] : seq N)).
635+ do 7 (case => [_|]; first by vm_cast_no_check (refl_equal ([::] : seq N))).
636+ by [].
637+ Time Qed .
638+
386639End Factorion.
387640
388- Check factorionE.
0 commit comments