Enigmatic Code

Programming Enigma Puzzles

Enigma 1536: List of primes

From New Scientist #2699, 14th March 2009 [link]

Use each of the digits 0 to 9 exactly twice to make a list of prime numbers. Your primes must all be different, and their sum must be as small as possible. Neither 0 nor 1 is a prime, and 0 may not be a leading digit.

What is the sum of the prime numbers in your list?

[enigma1536]

Advertisements

9 responses to “Enigma 1536: List of primes

  1. Jim Randell 11 June 2012 at 2:12 pm

    Here’s my original Perl solution – although you need to give it a maximum prime to search up to. If you already know the answer this is easy, 509 will find the minimal solutions. It takes 5m37s to do this. Searching up to the value of the minimum sum exhaustively verifies that there are no lower solutions. In this case checking primes up to 1287 takes 12m23s.

    use strict;
    
    # it turns out the solution is 1287, i.e.
    #  1287 = 2 5 23 47 61 67 83 89 401 509
    #  1287 = 2 5 29 47 61 67 83 89 401 503
    # so we can use that as the default for $MAX
    
    # although setting $MAX to 509 will find the same solutions faster
    # but isn't an exhaustive search (whereas 1287 is)
    
    my $MAX = $ARGV[0] || 509;
    
    my $MIN = [ 20 * $MAX ];
    
    # use a sieve to determine the non-primes
    
    my @NON_PRIME = (1, 1); # 0, 1 are non-prime
    
    my ($i, $j, $k, $v);
    my $t = time;
    
    for ($i = 2; $i <= $MAX; $i++) {
      next if $NON_PRIME[$i];
      for ($j = 2 * $i; $j <= $MAX; $j += $i) {
        $NON_PRIME[$j] ||= 1;
      }
    }
    
    # go through the list and set the values in the sieve to be the next highest prime
    $j = $MAX + 1;
    for ($i = $MAX; $i >=0; $i--) {
      if ($NON_PRIME[$i]) {
        $NON_PRIME[$i] = $j;
      } else {
        $j = $i;
      }
    }
    
    # now pick primes to use each digit twice
    check([(0) x 10], 2, 0);
    
    sub check {
      my ($d, $n, $sum, @primes) = @_;
    
      # see if we're done
      if ($d->[0] == 2 and $d->[1] == 2 and $d->[2] == 2 and $d->[3] == 2 and $d->[4] == 2 and
          $d->[5] == 2 and $d->[6] == 2 and $d->[7] == 2 and $d->[8] == 2 and $d->[9] == 2) {
        if ($sum <= $MIN->[0]) {
          $MIN = [ $sum, @primes ];
          print "$sum = @primes\n";
        }
        return;
      }
    
      # find the next candidate
      my %digits;
     TRY:
      # find the next prime
      $n = $NON_PRIME[$n] || $n;
      return if $sum + $n > $MIN->[0] or $n > $MAX;
    
      # check adding this in wouldn't use any digit more than twice
      %digits = ();
      for (split '', $n) {
        $digits{$_}++;
      }
      while (($k, $v) = each %digits) {
        if ($d->[$k] + $v > 2) {
          $n++;
          goto TRY;
        }
      }
    
      # OK add the candidate in
      while (($k, $v) = each %digits) {
        $d->[$k] += $v;
      }
    
      # and see if we can pick more primes to achieve a solution
      check($d, $n + 1, $sum + $n, @primes, $n);
    
      # then remove this prime
      while (($k, $v) = each %digits) {
        $d->[$k] -= $v;
      }
    
      # and try again from the next candidate
      check($d, $n + 1, $sum, @primes);
    }
    
    printf "MIN: %d = %s\n", $MIN->[0], join(' + ', @{$MIN}[1..@{$MIN}-1]);
    

    Solution: The sum of the prime numbers is 1287.

    • Jim Randell 11 June 2012 at 2:23 pm

      Here’s the same kind of approach in Python, using the Primes class from the enigma.py library. Again it tests up to a maximum value, for 509 it finds the solution in 1m36s (running under PyPy) and for 1287 it takes 3m22s.

      from enigma import Primes, split, printf
      
      # the smallest prime to have a zero in will be 101, so there must be
      # (at least) two primes in the list > 100.
      
      # allowable digits in a digit count
      valid_digits = set(('0', '1', '2'))
      
      # 509 will find the solution, 1287 is an exhaustive check
      import sys
      MAX = 509 if len(sys.argv) < 2 else eval(sys.argv[1])
      
      class Puzzle:
      
        def __init__(self, n):
          self.primes = list((p, sum(pow(10, i) for i in split(p, int))) for p in Primes(n))
          self.min = 20 * n
      
        def solve(self, dc=0, i=0, s=0, l=[]):
          # are we done?
          if dc == 2222222222:
            printf("sum={s} primes={l}")
            if s < self.min: self.min = s
            return
      
          # find the next candidate prime
          for j in range(i, len(self.primes)):
            p, d = self.primes[j]
            if s + p > self.min: return
            dcp = dc + d
            if not set(split(dcp)).issubset(valid_digits): continue
            # can we find a solution that includes this prime?
            self.solve(dcp, j + 1, s + p, l + [p])
          
      p = Puzzle(MAX)
      p.solve()
      printf("min = {p.min}")
      
      • Jim Randell 11 June 2012 at 2:24 pm

        And here’s a more comprehensive solution that increases the maximum value of primes it considers until it is greater than the minimum sum it has found, so it doesn’t need to know the answer before it starts. It runs in 1m11s (under PyPy).

        from enigma import Primes, split, printf
        
        # the smallest prime to have a zero in will be 101, so there must be
        # (at least) two primes in the list > 100.
        
        valid_digits = set(('0', '1', '2'))
        
        class Puzzle:
        
          def __init__(self, n=500, i=500):
            self.chunk = i
            self.sum = pow(10, 20)
            self.extend(n)
        
          def extend(self, n):
            self.max = n
            self.sieve = Primes(n)
            self.primes = self.sieve.list()[::-1]
            self.digits = list(sum(pow(10, i) for i in split(p, int)) for p in self.primes)
        
          def solve(self, ds=0, i=0, s=0, ps=[]):
            # are we done?
            if ds == 2222222222:
              if s < self.sum: self.sum = s
              printf("sum={s} primes={ps}", ps=list(reversed(ps)))
              return
            # find the next candidate prime
            for j in range(i, len(self.primes)):
              p = self.primes[j]
              # skip sums that are larger than the current smallest
              if s + p > self.sum: continue
              d = self.digits[j]
              # check we have not overused any of the digits
              if not set(split(ds + d)).issubset(valid_digits): continue
              self.solve(ds + d, j + 1, s + p, ps + [p])
        
          def run(self):
            # the largest prime must be 107 or more
            l = len(self.primes) - self.primes.index(107)
            while True:
              # choose the largest prime
              for i in range(l, len(self.primes)+1):
                p = self.primes[-i]
                # if the prime is larger than current sum we've completed the exhaustive search
                if not(p < self.sum):
                  printf("MIN = {self.sum}")
                  return
                self.solve(self.digits[-i], len(self.primes) - i + 1, p, [p])
              # if we've not reached the end of the search extend the primes list
              n = self.max + self.chunk
              l = len(self.primes) + 1
              printf("extending max prime from {self.max} to {n}")
              self.extend(n)
        
        
        p = Puzzle()
        p.run()
        
  2. arthurvause 27 April 2013 at 10:11 am

    For some reason, I didn’t produce a solution for this one at the time it appeared, but it is quite an interesting problem.

    Here is a recursive solution, running in 6 seconds on my laptop with Python 2.7

    from  itertools import combinations as combs
    
    def notexceed(s):
      set_s = set(s)
      if " " in set_s:
        set_s.remove(" ")
      for n in set_s:
          if s.count(n) > 2:
              return False
      return True
    
    primes = set([2]+range(3,1000,2))-set.union(*[set(range(n*n,1000,n*2)) for n in xrange(3,32,2)])
    
    primescontainingsingle = [0]*10
    primescontainingdouble = [0]*10
    for n in xrange(10):
        primescontainingsingle[n] = { str(x) for x in primes if str(x).count(str(n))==1 }
        primescontainingdouble[n] = { str(x) for x in primes if str(x).count(str(n))==2 }
    
    seq = [0,8,2,4,5,6,9,3,1,7]
    
    def add_digit(n,s,total):
      global results
      global besttotal
      
      if total>besttotal:
        return
    
      L = len(s)-s.count(" ")
      if L>20:
          return
    
      if L==20:
        if notexceed(s):
          besttotal=total
          print "result", s, ", Total=",total
        return
    
      if n==10:
        return  
    
      if notexceed(s):
    
        d=seq[n]
        d_sofar = s.count(str(d))
        if d_sofar==0:
          for p,q in combs(primescontainingsingle[d],2):
            if not(" "+p+" " in s) and not(" "+q+" " in s): 
              add_digit(n+1,s+p+" "+q+" ",total+int(p)+int(q))
    
          for p in primescontainingdouble[d]:
             if not(" "+p+" " in s): 
              add_digit(n+1,s+p+" ",total+int(p))
                  
        elif d_sofar==1:  
          for p in primescontainingsingle[d]:
             if not(" "+p+" " in s): 
              add_digit(n+1,s+p+" ",total+int(p))
    
        elif d_sofar==2:  
          add_digit(n+1,s,total)
    
        else:
          return
    
    besttotal=100000
    
    add_digit(0,' ',0)
    
    
  3. Jim Randell 27 April 2013 at 12:50 pm

    I had another look at this one, and recoded my previous best algorithm to bring the GLPK solver to bear on it (using PyMathProg, which I discovered when writing a solution for Enigma 82). This program runs in 495ms (although if you start it with a max prime of 1287 this goes down to 229ms, solve(509) by itself only takes 112ms).

    import pymprog
    from enigma import Primes, irange, printf
    
    # digits
    digits = list(irange(0, 9))
    
    # find a solution for primes up to maxp
    def solve(maxp):
    
      primes = list(Primes(maxp))
    
      # map primes to digit counts
      dc = dict()
      for p in primes:
        s = str(p)
        dc[p] = list(s.count(str(d)) for d in digits)
    
      # create the model
      m = pymprog.model('enigma1536')
    
      # we need a descision variable for each prime
      x = m.var(primes, 'x', bool)
    
      # objective: minimise the sum of the primes used
      m.min(sum(p * x[p] for p in primes))
    
      # constraints: each digit is used twice
      m.st([sum(dc[p][d] * x[p] for p in primes) == 2 for d in digits])
    
      # find a solution
      m.solve(int)
    
      # which primes are used?
      ps = list(p for p in primes if bool(x[p].primal))
      s = sum(ps)
    
      printf("sum={s} {ps}")
      return s
    
    import sys
    n = 500 if len(sys.argv) < 2 else int(sys.argv[1])
    i = 500 if len(sys.argv) < 3 else int(sys.argv[2])
    
    # increase the max prime until it's larger than the best solution
    while True:
      s = solve(n)
      if not(n < s): break
      n += i
    
    printf("min sum = {s}")
    
  4. Paul 13 April 2015 at 3:19 pm

    Here is a slightly different method using MMa. Experimentation has shown that any other than 2 x 1 digit and 6 x 2 digit primes doesn’t produce any solutions, so using that information the program generates all subsets of separate 1 and 2 digit primes which gives a list of 325584 possible solutions. from this list it is reduced to 2681 possible solutions that would leave 6 digits from the possible remaining 20 digits. These 6 digits are permuted and separated in to 2 x 3 digit numbers and both tested to be prime, if so we have a solution. The program finds all 26 solutions with an upper prime limit of 997, you can choose from the list the 2 minimum solutions.

    Timing[a=Prime[Range[4]];
    b=Prime[Range[5,25]];
    d=Subsets[a,{2}];
    e=Subsets[b,{6}];
    gg=Tuples[{d,e}];hh=Cases[gg,{{a_,b_},{d_,e_,f_,g_,h_,i_}}/;
    Max[Last/@Tally[Flatten[IntegerDigits[Flatten[{a,b,d,e,f,g,h,i}]]]]]==2];
    Do[ww=RotateRight[DigitCount[FromDigits[Flatten[IntegerDigits/@Flatten[hh[[qw]]]]]]];lis={};
    Do[If[ww[[q]]<2,ee=2-ww[[q]];Do[AppendTo[lis,q-1],{gt,1,ee}]],{q,1,10}];
    u=Cases[Permutations[lis],{a_,b_,c_,d_,e_,f_}/;
    a!=0&&d!=0&&c!=0&&f!=0&&PrimeQ[t=100a+10b+c]&&PrimeQ[k=100d+10e+f]&&t<k];
    u=Flatten[u];
    If[Length[u]>0,Print[hh[[qw]],"  ",FromDigits[u[[1;;3]]],"  ",FromDigits[u[[4;;6]]],"  ",FromDigits[u[[1;;3]]]+FromDigits[u[[4;;6]]]+Total[Total/@hh[[qw]]]]],{qw,1,Length[hh]}]]
    

    it runs in 7.1 seconds.

Leave a Comment

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: