(in-package "ACL2")
;; This book illustrates a simple application of the theory of quadratic
;; residues to the Mersenne prime problem. We show that certain Mersenne
;; numbers may be factored by applying a theorem of Euler, which combines
;; two results: (1) Euler's criterion for quadratic residues and (2) the
;; second supplement to the law of quadratic reciprocity, which computes
;; the quadratic character of 2 modulo an odd prime p.
;; The two theorems and the relevant definitions are contained in the
;; book "gauss":
(include-book "gauss")
(comp t)
; We skip some forms by default, because they take a long time to run.
(defmacro maybe-skip (x)
; Redefine to be x if you want to execute forms inside (maybe-skip ...).
(declare (ignore x))
'(value-triple nil))
;; For a prime p, the Mersenne number Mp is defined to be 2^p-1:
(defun m (p) (1- (expt 2 p)))
;; We would like to determine the primality of Mp. For small values of p,
;; this can be done simply by executing the function primep:
(defthm mersenne-19
(primep (m 19)))
;;[Time: .4 seconds]
(maybe-skip
(defthm mersenne-31
(primep (m 31)))
)
;;[Time: 65 minutes]
;; The timing of mersenne-31 suggests that the computation of mersenne-61
;; is intractable:
;;(defthm mersenne-61
;; (primep (m 61)))
;;[Time: a billion hours?]
;; Naturally, we can handle larger composites than primes:
(defthm mersenne-23
(not (primep (m 23))))
;;[Time: .02 seconds]
(maybe-skip
(defthm mersenne-67
(not (primep (m 67))))
)
;;[Time: 288 seconds]
(maybe-skip
(defthm mersenne-999671
(not (primep (m 999671))))
)
;;[Time: 165 minutes]
;;(defthm mersenne-19876271
;; (not (primep (m 19876271))))
;;[Error: Attempt to create an integer that is too large to represent.]
;; There is an obvious optimization of this procedure, based on
;; the observation that if n has a proper divisor, then it has one
;; that is at most sqrt(n). Thus, we define an alternative to
;; the function least-divisor that stops at sqrt(n), and
;; establish a rewrite rule:
(defun least-divisor-2 (k n)
(declare (xargs :measure (nfix (- n k))))
(if (and (integerp n)
(integerp k)
(< 1 k)
(<= k n))
(if (> (* k k) n)
n
(if (divides k n)
k
(least-divisor-2 (1+ k) n)))
nil))
(comp t)
(defthm least-divisor-rewrite-lemma-1
(implies (and (integerp n)
(integerp k)
(< 1 k)
(<= k n)
(<= (* (least-divisor k n) (least-divisor k n)) n))
(equal (least-divisor k n)
(least-divisor-2 k n)))
:rule-classes ()
:hints (("Subgoal *1/7" :use (least-divisor-divides))
("Subgoal *1/1" :use (least-divisor-divides))))
(defthm least-divisor-rewrite-lemma-2
(implies (and (integerp n)
(> n 1)
(not (= n (least-divisor 2 n))))
(<= (* (least-divisor 2 n) (least-divisor 2 n)) n))
:rule-classes ()
:hints (("Goal" :in-theory (enable divides)
:use ((:instance least-divisor-divides (k 2))
(:instance least-divisor-is-least (k 2) (d (/ n (least-divisor 2 n))))
(:instance *-weakly-monotonic (x (least-divisor 2 n)) (y 1) (y+ (/ n (least-divisor 2 n))))
(:instance *-weakly-monotonic (x (least-divisor 2 n)) (y+ 1) (y (/ n (least-divisor 2 n))))))))
(defthm least-divisor-rewrite-lemma-3
(implies (and (integerp n)
(integerp k)
(< 1 k)
(<= k n)
(= n (least-divisor k n)))
(= n (least-divisor-2 k n)))
:rule-classes ()
:hints (("Goal" :use (least-divisor-divides
(:instance least-divisor-is-least (d (least-divisor-2 2 n)))))))
(defthm least-divisor-rewrite
(equal (least-divisor 2 n)
(least-divisor-2 2 n))
:hints (("Goal" :use (least-divisor-rewrite-lemma-2
(:instance least-divisor-rewrite-lemma-3 (k 2))
(:instance least-divisor-rewrite-lemma-1 (k 2))))))
;; In order to arrange that the more efficient variant is executed, we must
;; disable the executable counterparts of primep and least-proper-divisor
;; and enable the logical function primep. This allows us to establish
;; primality somewhat more efficiently:
(in-theory (enable primep))
(in-theory (disable (primep) (least-divisor)))
(defthm mersenne-31-revisited
(primep (m 31)))
;;[Time: .05 seconds]
(maybe-skip
(defthm mersenne-61
(primep (m 61)))
)
;;[Time: 54 minutes]
;; However, the above optimization does not help at all in the composite case.
;; In the special case where mod(p,4) = 3 and 2*p+1 is prime, the following
;; theorem may be used. This theorem is attributed to Euler in Hardy and
;; Wright, where it appears as Theorem 103. It is an immediate comsequence
;; of Euler's Criterion and the Second Supplement.
(defthm theorem-103
(implies (and (primep p)
(= (mod p 4) 3)
(primep (1+ (* 2 p))))
(divides (1+ (* 2 p)) (m p)))
:rule-classes ()
:hints (("Goal" :use ((:instance second-supplement (p (1+ (* 2 p))))
(:instance euler-criterion (m 2) (p (1+ (* 2 p))))
(:instance mod-sum (a 1) (b (* 2 p)) (n 8))
(:instance mod-prod (k 2) (m p) (n 4))
(:instance divides-leq (x (1+ (* 2 p))) (y 2))
(:instance divides-mod-equal (a (expt 2 p)) (b 1) (n (1+ (* 2 p))))))))
;; With the exception of p = 3 (in which case Mp = 2*p+1), Mp is composite
;; for any p satisfying the hypothesis of theorem-103:
(defthm expt-2-p-bnd
(implies (and (integerp p)
(> p 3))
(> (expt 2 p) (* 2 (1+ p))))
:rule-classes ())
(defthm theorem-103-corollary
(implies (and (primep p)
(= (mod p 4) 3)
(> p 3)
(primep (1+ (* 2 p))))
(not (primep (m p))))
:hints (("Goal" :use (theorem-103
expt-2-p-bnd
(:instance primep-no-divisor (p (1- (expt 2 p))) (d (1+ (* 2 p))))))))
;; In order to arrange for theorem-103-corollary to be used as a
;; rewrite rule, we must disable the logical expansion of both
;; primep and m. We must also re-enable the executable counterpart
;; of primep to allow the computation of primep(p) and primep(2*p+1):
(in-theory (disable primep m (m)))
(in-theory (enable (primep)))
(defthm mersenne-23-revisited
(not (primep (m 23))))
;;[Time: .01 seconds]
(defthm mersenne-999671-revisited
(not (primep (m 999671))))
;;[Time: .01 seconds]
(maybe-skip
(defthm mersenne-19876271
(not (primep (m 19876271))))
)
;;[Time: 47 seconds]