guile-project-euler-solutions/027/solution.scm

;; https://projecteuler.net/problem=27

;; Quadratic primes
;; Problem 27

;; Euler discovered the remarkable quadratic formula:

;; n^2 + n + 41

;; It turns out that the formula will produce 40 primes for
;; the consecutive integer values

;; 0 <= n <= 39

;; . However, when

;; n = 40
;; 40^2 + 40 + 41 = 40(40 + 1) + 41

;; is divisible by 41, and certainly when

;; n = 41
;; 41^2 + 41 + 41

;; is clearly divisible by 41.

;; The incredible formula

;; n^2 - 79n + 1601

;; was discovered, which produces 80 primes for the
;; consecutive values

;; 0 <= n <= 79

;; . The product of the coefficients, -79 and 1601, is
;; -126479.

;; Considering quadratics of the form:

;;     n^2 + an + b, where |a| < 1000 and |b| <= 1000

;;     where |n| is the modulus/absolute value of n
;;     e.g. |11| = 11 and |-4| = 4

;; Find the product of the coefficients, a and b, for the
;; quadratic expression that produces the maximum number of
;; primes for consecutive values of n, starting with n = 0.


(import
 (except (rnrs base) let-values map)
 (only (guile)
       lambda* λ)
 (ice-9 match)
 (contract)
 (prefix (lib math) math:)
 (lib naive-prime-test)
 (lib print-utils)
 (lib parallelism)
 (lib segment))


(define make-quadratic-formula
  (λ (coeff-a coeff-b)
    ;; (when (and (= coeff-a -61) (= coeff-b 971))
    ;;   (print "making quadratic formula with a =" coeff-a "and b =" coeff-b))
    (λ (num)
      (+ (math:square num)
         (* coeff-a num)
         coeff-b))))


;; (define-with-contract euler-remarkable-formula
;;   (require (integer? num)
;;            (and (>= num 0) (< num 40)))
;;   (ensure (integer? <?>)
;;           (positive? <?>)
;;           (prime? <?>))
;;   (λ (num)
;;     ((make-quadratic-formula 1 41) num)))


;; (define-with-contract incredible-formula
;;   (require (integer? num)
;;            (and (>= num 0) (< num 80)))
;;   (ensure (integer? <?>)
;;           (positive? <?>)
;;           (prime? <?>))
;;   (λ (num)
;;     ((make-quadratic-formula -79 1601) num)))


(define-with-contract number-of-primes
  (require (procedure? formula))
  (ensure (math:natural-number? <?>))
  (λ (formula)
    (let iter ([num 0])
      ;; (print "trying formula for num:" num "which is:" (formula num))
      (cond
       [(prime? (formula num)) (iter (+ num 1))]
       [else num]))))


;; (number-of-primes (make-quadratic-formula -61 971)) is overlooked?

(define find-most-primes-coeffs
  (λ (a-start a-end)
    (print "starting search for" a-start "<= a <=" a-end)
    (let ([limit-abs-b 1000])
      (let iter ([coeff-a a-start]
                 [coeff-b (* -1 limit-abs-b)]
                 [maximum-num-primes 0]
                 [result (cons 0 (cons 0 0))])
        ;; (when (and (= coeff-a -61) (= coeff-b 971))
        ;;   (print "a =" coeff-a
        ;;          "and b =" coeff-b))
        (cond
         [(> coeff-b limit-abs-b)
          ;; increase a, reset b
          ;; (print "increment a to" (+ coeff-a 1))
          (iter (+ coeff-a 1)
                (* -1 limit-abs-b)
                maximum-num-primes
                result)]

         [(> coeff-a a-end)
          ;; (print "reached end of segment [" a-start ";" a-end "]" ", result is:" result)
          result]

         [else
          (let* ([formula (make-quadratic-formula coeff-a coeff-b)]
                 [num-primes (number-of-primes formula)])
            (cond
             [(> num-primes maximum-num-primes)
              ;; (print "new maximum of produced primes for"
              ;;        "a:" coeff-a
              ;;        "b:" coeff-b
              ;;        "->" num-primes)
              (iter coeff-a
                    (+ coeff-b 1)
                    num-primes
                    (cons num-primes (cons coeff-a coeff-b)))]

             [else
              (iter coeff-a
                    (+ coeff-b 1)
                    maximum-num-primes
                    result)]))])))))


(define primes-and-coeffs
  (let ([limit-abs-a 999])
    (let* ([num-cores 32]
           [segments (segment (* -1 limit-abs-a)
                              limit-abs-a
                              num-cores)])
      ;; (print "segments:" segments)
      (run-in-parallel segments
                       (λ (seg)
                         (find-most-primes-coeffs (segment-start seg)
                                                  (segment-end seg)))
                       (λ (acc current)
                         ;; (print "comparing result" acc "with result" current)
                         (cond
                          [(> (car acc) (car current))
                           acc]
                          [else current]))
                       -inf.0))))


(match primes-and-coeffs
  [(num-of-primes a . b)
   (print "quadratic formula with coefficients"
          "a =" a "and"
          "b =" b
          "has" num-of-primes "primes")
   (print "result:" (* a b))])