大学的时候挺喜欢解 Project Euler 上的题目的，尤其是它不在乎答题者使用哪一门编程语言，甚至还有很多参与者是使用 pen&paper 来解题的。去年开始重新开始做 Project Euler 上的题目，而第 69 题则是最近刚刚解决的一题。惭愧的是，因为不晓得欧拉函数的计算公式（甚至都没有想过欧拉函数有没有可以用来计算的公式），所以这一题我是用暴力计算的方法来解决的。尽管花了 40 分钟左右才找出了问题的答案，但欧拉函数的计算方法本身还是让我觉得挺有意思的，下面我就来讲讲我在计算欧拉函数方面做的一些尝试。
69 题本身很容易读懂，就是要找到一个不大于一百万的正整数 n，这个 n 与以它为参的欧拉函数值的比值能够达到最大。欧拉函数的介绍可以看维基百科，简而言之，就是在不大于 n 的正整数中与 n 互质的数的个数，具体的例子可以参见 69 题描述中给出的表格。
网上可以找到别人的解法，基本的思路是：按从小到大的顺序，对于一百万以内的每个素数，都计算出它们的倍数的欧拉函数值的一部分——即对于素数 p 计算出 1 -1/ p 并与这个位置上原来的值相乘。当遍历完一百万以内的所有素数后，也就计算出每一个位置上的欧拉函数值，再遍历一次就可以计算出比值最大的数字了。
但我今天要讲的是笨方法。
笨方法的关键就是乖巧地计算出每一个数字的欧拉函数值。而其中最笨的，当属挑出每一个不大于 n 的因数，计算它们与 n 的最大公约数，并根据这个最大公约数是否为 1，来决定是否给某一个计数器加一。示例代码如下
(defun phi (n)
(let ((count 0))
(dotimes (i n count)
(when (= (gcd (1+ i) n) 1)
(incf count)))))
这个 phi 可以稍微改进一下。例如，如果一个数 a 与 n 不是互质的，那么 a 的倍数（小于 n 的那一些）也一定不会与 n 互质。因此，每当遇到这么一个因数 a，就知道后续的 2a、3a 等等，都不再需要计算其与 n 的最大公约数了。改进后的算法代码如下
(defun phi (n)
“ 通过将不互质的比特设置为 1 并计算为 0 的比特的个数来计算 phi 函数 ”
(let ((bits (make-array n :element-type ‘bit :initial-element 0))
(count 0))
(dotimes (i n)
(cond ((= (bit bits i) 1)
;; 该比特已经为 1，说明已经在比它小的倍数被处理时一并被标记了
)
((= i (1- n))
;; 只处理比上界要小的数字
)
((/= (gcd (1+ i) n) 1)
;; 除了当前这个不互质的数字之外，还需要将这个数字的倍数也一并处理
(dotimes (j (floor n (1+ i)))
(let* ((j (1+ j))
(m (* (1+ i) j)))
(setf (bit bits (1- m)) 1))))
(t (incf count))))
count))
为了节省内存空间，这里用了一个 bitmap 来标记小于 n 的每一个数字是否与 n 互质——1 表示不互质，0 表示互质。
其实并不需要遍历所有比 n 小的数字，只要遍历所有 n 的素因子即可。比如，将 8 分解为素因子，就是 3 个 2，那么对于小于 8 的所有 2 的倍数（4 和 6）都是与 8 不互质的。基于这个方法，将所有的素因子的倍数所对应的位置为 1，再数一下总共有多少个 0 比特即可。
对每个 n 都进行质因数分解效率不高，先生成一个一百万以内的素数表吧
(defun primep (n)
(cond ((< n 2) nil)
((= 2 n) t)
(t
(let ((bnd (truncate (sqrt n))))
(labels ((rec (test)
(cond ((> test bnd) t)
((= 0 (rem n test)) nil)
(t (rec (1+ test))))))
(rec 2))))))
(defvar *primes-in-1000000* nil)
(defun generate-primes-in-1000000 ()
(dotimes (i 1000000)
(when (primep (1+ i))
(push (1+ i) *primes-in-1000000*)))
(setf *primes-in-1000000* (nreverse *primes-in-1000000*)))
然后对于一个给定的 n，遍历所有小于它的素数并对相应的倍数所在的比特置一就可以了，示例代码如下
(defun phi3 (n)
“ 直接用素数表来做筛法 ”
(prog
((bits (make-array n :element-type ‘bit :initial-element 0)))
(dolist (num *primes-in-1000000*)
(cond ((> num n)
(go :return))
((zerop (mod n num))
(labels ((aux (i)
(when (< i n)
(setf (bit bits i) 1)
(aux (+ i num)))))
(aux (1- num))))))
:return
(return-from phi3
(count-if (lambda (bit) (zerop bit)) bits))))
PS：写这个 phi3 的时候发现 Common Lisp 提供了一个 prog 宏，这个宏倒是真的挺好用。
改进了两轮，其实这仍然是笨方法。即便是用 phi3，用来计算题目的答案也花了 40 多分钟。
全文完