Euclid-Mullin sequence

From Rosetta Code
Revision as of 13:38, 28 October 2021 by PureFox (talk | contribs) (Created a new draft task, Euclid-Mullin sequence, and added a Wren solution.)
(diff) ← Older revision | Latest revision (diff) | Newer revision → (diff)
Euclid-Mullin sequence is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.
Definition

The Euclid–Mullin sequence is an infinite sequence of distinct prime numbers, in which each element is the least prime factor of one plus the product of all earlier elements.

The first element is usually assumed to be 2. So the second element is : (2) + 1 = 3 and the third element is : (2 x 3) + 1 = 7 as this is prime.

Although intermingled with smaller elements, the sequence can produce very large elements quite quickly and only the first 51 have been computed at the time of writing.

Task

Compute and show here the first 16 elements of the sequence of, if your language does not support arbitrary precision arithmetic, as many as you can.

Stretch goal

Compute the next 11 elements of the sequence.

Reference

OEIS sequence A000945

Wren

This uses the Pollard Rho algorithm to try and speed up the factorization of the 15th element but overall time still slow at around 32 seconds. <lang ecmascript>import "./big" for BigInt

var zero = BigInt.zero var one = BigInt.one var two = BigInt.two var ten = BigInt.ten var max = BigInt.new(100000)

var pollardRho = Fn.new { |n, c|

   var g = Fn.new { |x, y| (x*x + c) % n }
   var x = two
   var y = two
   var z = one
   var d = max + one
   var count = 0
   while (true) {
       x = g.call(x, n)
       y = g.call(g.call(y, n), n)
       d = (x - y).abs % n
       z = z * d
       count = count + 1
       if (count == 100) {
           d = BigInt.gcd(z, n)
           if (d != one) break
           z = one
           count = 0
       }
   }
   if (d == n) return zero
   return d

}

var smallestPrimeFactorWheel = Fn.new { |n|

   if (n.isProbablePrime(5)) return n
   if (n % 2 == zero) return BigInt.two
   if (n % 3 == zero) return BigInt.three
   if (n % 5 == zero) return BigInt.five
   var k = BigInt.new(7)
   var i = 0
   var inc = [4, 2, 4, 2, 4, 6, 2, 6]
   while (k * k <= n) {
       if (n % k == zero) return k
       k = k + inc[i]
       if (k > max) return null
       i = (i + 1) % 8
   }

}

var smallestPrimeFactor = Fn.new { |n|

   var s = smallestPrimeFactorWheel.call(n)
   if (s) return s
   var c = one
   s = n
   while (n > max) {
       var d = pollardRho.call(n, c)
       if (d == 0) {
           if (c == ten) Fiber.abort("Pollard Rho doesn't appear to be working.")
           c = c + one           
       } else {
           // can't be sure PR will find the smallest prime factor first
           s = BigInt.min(s, d)
           n = n / d
           if (n.isProbablePrime(2)) return BigInt.min(s, n)
       }
   }
   return s

}

var k = 16 System.print("First %(k) terms of the Euclid–Mullin sequence:") System.print(2) var prod = BigInt.two var count = 1 while (count < k) {

   var t = smallestPrimeFactor.call(prod + one)
   System.print(t)
   prod = prod * t
   count = count + 1

} </lang>

Output:
First 16 terms of the Euclid–Mullin sequence:
2
3
7
43
13
53
5
6221671
38709183810571
139
2801
11
17
5471
52662739
23003