# 5e5dd3fd5ae0

Problems 71 and 73

author | Steve Losh <steve@stevelosh.com> |
---|---|

date | Thu, 01 Nov 2018 19:44:26 -0400 |

parents | 99fa06464617 |

children | 2a837c783fe7 |

branches/tags | (none) |

files | src/problems.lisp src/utils.lisp |

## Changes

--- a/src/problems.lisp Tue Jul 31 16:18:33 2018 +0000 +++ b/src/problems.lisp Thu Nov 01 19:44:26 2018 -0400 @@ -1592,6 +1592,60 @@ (while (<= n 1000000)) (finally (return result)))) +(defun problem-71 () + ;; Consider the fraction, n/d, where n and d are positive integers. If n<d and + ;; HCF(n,d)=1, it is called a reduced proper fraction. + ;; + ;; If we list the set of reduced proper fractions for d ≤ 8 in ascending order + ;; of size, we get: + ;; + ;; 1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, + ;; 5/7, 3/4, 4/5, 5/6, 6/7, 7/8 + ;; + ;; It can be seen that 2/5 is the fraction immediately to the left of 3/7. + ;; + ;; By listing the set of reduced proper fractions for d ≤ 1,000,000 in + ;; ascending order of size, find the numerator of the fraction immediately to + ;; the left of 3/7. + (iterate + ;; We could theoretically iterate through F_1000000 til we find 3/7, but + ;; we'd have about 130000000000 elements to churn though. We can do it much + ;; faster. + ;; + ;; Instead, we start by computing F₇ and finding 3/7 and whatever's to the + ;; left of it. We can find what's in between these two in the next + ;; iteration by taking their mediant, then just repeat until we hit the + ;; limit. + (with limit = 1000000) + (with target = 3/7) + (with guess = (iterate (for x :in-farey-sequence (denominator target)) + (for px :previous x) + (finding px :such-that (= x target)))) + (for next = (mediant guess target)) + (if (> (denominator next) limit) + (return (numerator guess)) + (setf guess next)))) + +(defun problem-73 () + ;; Consider the fraction, n/d, where n and d are positive integers. If n<d and + ;; HCF(n,d)=1, it is called a reduced proper fraction. + ;; + ;; If we list the set of reduced proper fractions for d ≤ 8 in ascending order + ;; of size, we get: + ;; + ;; 1/8, 1/7, 1/6, 1/5, 1/4, 2/7, 1/3, 3/8, 2/5, 3/7, 1/2, 4/7, 3/5, 5/8, 2/3, + ;; 5/7, 3/4, 4/5, 5/6, 6/7, 7/8 + ;; + ;; It can be seen that there are 3 fractions between 1/3 and 1/2. + ;; + ;; How many fractions lie between 1/3 and 1/2 in the sorted set of reduced + ;; proper fractions for d ≤ 12,000? + (iterate + ;; Just brute force this one with our fancy Farey sequence iterator. + (for x :in-farey-sequence 12000) + (while (< x 1/2)) + (counting (< 1/3 x)))) + (defun problem-74 () ;; The number 145 is well known for the property that the sum of the factorial ;; of its digits is equal to 145: @@ -2365,6 +2419,8 @@ (test p63 (is (= 49 (problem-63)))) (test p67 (is (= 7273 (problem-67)))) (test p69 (is (= 510510 (problem-69)))) +(test p71 (is (= 428570 (problem-71)))) +(test p73 (is (= 7295372 (problem-73)))) (test p74 (is (= 402 (problem-74)))) (test p79 (is (= 73162890 (problem-79)))) (test p81 (is (= 427337 (problem-81))))

--- a/src/utils.lisp Tue Jul 31 16:18:33 2018 +0000 +++ b/src/utils.lisp Thu Nov 01 19:44:26 2018 -0400 @@ -167,6 +167,30 @@ (cl-strings:split l ,delimiter%)))) (,kwd ,var :next (parse-line (next-line)))))))) +(defmacro-driver (FOR var IN-FAREY-SEQUENCE n) + "Iterate `var` through the Farey sequence Fₙ in ascending order. + + The Farey sequence Fₙ is the sorted sequence of all reduced fractions between + 0 and 1 (inclusive) with denominators ≤ `n`. + + Example: + + (iterate (for x :in-farey-sequence 5) (collect x)) + ; => + (0 1/5 1/4 1/3 2/5 1/2 3/5 2/3 3/4 4/5 1) + + " + (with-gensyms (n% prev) + (let ((kwd (if generate 'generate 'for))) + `(progn + (with ,n% = ,n) + (with ,prev = 0) + (,kwd ,var :do-next + (if-first-time + (setf ,var 0) + (psetf ,var (or (next-farey-term ,prev ,var ,n%) (terminate)) + ,prev ,var))))))) + (defmacro when-first-time (&body body) `(if-first-time @@ -345,6 +369,16 @@ (u (- 1 v w))) (values u v w))) +(defun mediant (x y) + "Return the mediant of `x` and `y`." + ;; The mediant operation is: + ;; + ;; A C A + C + ;; - ⊕ - = ----- + ;; B D B + D + (/ (+ (numerator x) (numerator y)) + (+ (denominator x) (denominator y)))) + ;;;; Digits ------------------------------------------------------------------- (defun digits-length (n &optional (radix 10)) @@ -781,8 +815,6 @@ (define-modify-macro adjoinf (item &rest keyword-args) adjoin%) - - (defun pythagorean-triplet-p (a b c) (math a ^ 2 + b ^ 2 = c ^ 2)) @@ -818,6 +850,40 @@ (recur (mv* d triple)))))) +;;;; Farey Sequences ---------------------------------------------------------- +(defun next-farey-term (x y n) + "Return the next term in the Farey sequence Fₙ after `x` and `y`, or `nil`. + + For convenience, if `y` is zero `x` is not examined. + + If `y` is `1` the sequence is complete and `nil` is returned. + + " + (declare (type (rational 0 1) x y) + (type (integer 1) n)) + (case y + (0 (/ 1 n)) + (1 nil) + (t (let* ((a (numerator x)) + (b (denominator x)) + (c (numerator y)) + (d (denominator y)) + (k (floor (+ n b) d)) + (p (- (* k c) a)) + (q (- (* k d) b))) + (/ p q))))) + +(defun farey (n) + "Return the Farey sequence Fₙ as a fresh list." + (iterate (for x :in-farey-sequence n) + (collect x))) + +(defun approximate-farey-length (n) + "Return an approximation of the length of the Farey sequence Fₙ." + (/ (* 3 (expt n 2)) + (expt pi 2))) + + ;;;; Geometric Numbers -------------------------------------------------------- (defun squarep (n) "Return whether `n` is a perfect square." @@ -974,7 +1040,6 @@ (* precision (round number precision))) - ;;;; A* Search ---------------------------------------------------------------- (defstruct path state