extra/project-euler/027/027.factor

   1 ! Copyright (c) 2008 Aaron Schaefer.
   2 ! See https://factorcode.org/license.txt for BSD license.
   3 USING: kernel math math.primes project-euler.common sequences ;
   4 IN: project-euler.027
   5
   6 ! https://projecteuler.net/problem=27
   7
   8 ! DESCRIPTION
   9 ! -----------
  10
  11 ! Euler published the remarkable quadratic formula:
  12
  13 !     n² + n + 41
  14
  15 ! It turns out that the formula will produce 40 primes for the
  16 ! consecutive values n = 0 to 39. However, when n = 40, 402 + 40
  17 ! + 41 = 40(40 + 1) + 41 is divisible by 41, and certainly when
  18 ! n = 41, 41² + 41 + 41 is clearly divisible by 41.
  19
  20 ! Using computers, the incredible formula n² - 79n + 1601 was
  21 ! discovered, which produces 80 primes for the consecutive
  22 ! values n = 0 to 79. The product of the coefficients, -79 and
  23 ! 1601, is -126479.
  24
  25 ! Considering quadratics of the form:
  26
  27 !     n² + an + b, where |a| < 1000 and |b| < 1000
  28
  29 !     where |n| is the modulus/absolute value of n
  30 !     e.g. |11| = 11 and |-4| = 4
  31
  32 ! Find the product of the coefficients, a and b, for the
  33 ! quadratic expression that produces the maximum number of
  34 ! primes for consecutive values of n, starting with n = 0.
  35
  36
  37 ! SOLUTION
  38 ! --------
  39
  40 ! b must be prime since n = 0 must return a prime
  41 ! a + b + 1 must be prime since n = 1 must return a prime
  42 ! 1 - a + b must be prime as well, hence >= 2. Therefore:
  43 !    1 - a + b >= 2
  44 !        b - a >= 1
  45 !            a < b
  46
  47 <PRIVATE
  48
  49 : source-027 ( -- seq )
  50     1000 <iota> [ prime? ] filter [ dup [ neg ] map append ] keep
  51     cartesian-product concat [ first2 < ] filter ;
  52
  53 : quadratic ( b a n -- m )
  54     dup sq -rot * + + ;
  55
  56 : (consecutive-primes) ( b a n -- m )
  57     3dup quadratic prime? [ 1 + (consecutive-primes) ] [ 2nip ] if ;
  58
  59 : consecutive-primes ( a b -- m )
  60     swap 0 (consecutive-primes) ;
  61
  62 : max-consecutive ( seq -- elt n )
  63     dup [ first2 consecutive-primes ] map dup supremum
  64     over index [ swap nth ] curry bi@ ;
  65
  66 PRIVATE>
  67
  68 : euler027 ( -- answer )
  69     source-027 max-consecutive drop product ;
  70
  71 ! [ euler027 ] 100 ave-time
  72 ! 111 ms ave run time - 6.07 SD (100 trials)
  73
  74 ! TODO: generalize max-consecutive/max-product (from #26) into a new word
  75
  76 SOLUTION: euler027