Yesterday’s post from Programming Praxis poses an interesting problem: find the largest prime *n* such that the result of repeatedly removing each digit of *n* from left to right is also always prime.

For example, 6317 would be such a number, as not only is it prime, but so are 317, 17, and 7.

This is actually a surprisingly straightforward problem, if instead of starting at the largest such number and checking to make sure that each part along the way is prime, you start with the single digit primes (2, 3, 5, and 7) and repeatedly add on numbers to the left so long as they are prime.

For example, 3 would become 13, 23, 43, 53, 73, and 83. 13 in turn would become 113, 313, and 613. It turns out though that this list doesn’t grow unbounded. Eventually each branch will terminate (a proof for why this happens would be interesting, but I’m not sure how it would work–or even if such a thing is provable short of proof by enumerating all of the solutions).

Enough on discussion though, it’s code time!

(Only Racket today as I’m running a bit short on time for my writeup. If anyone would like a Python version as well, let me know in the comments and I’ll write one up.)

; find the largest prime such that each number produced by ; removing digits from the left side is still prime ; #:root - the prime to start at, 0 for all primes ; #:limit - the largest prime to check (define (largest-nested-prime #:root [root 0] #:limit [limit +inf.0]) ; loop starting at the root (let loop ([n root]) (cond ; if we're over the limit, no largest prime [(> n limit) 0] ; otherwise, try each prefixed digit while still prime [else (define multiplier (expt 10 (digits n))) ; use for/fold to emulate what for/max would do (for/fold ([best n]) ([i (in-range 1 10)] #:when (prime? (+ (* multiplier i) n))) ; return the best of the current and the nested (max best (loop (+ (* multiplier i) n))))])))

That’s a whopping 10 lines of code if you remove the comments and collapse the define back into a single line. Not too bad for what it’s calculating.

It takes rather a while to run, but eventually you’ll get an answer:

> (largest-nested-prime) 357686312646216567629137

What I’m curious about though, is there a better way to do this? I used the trial division method of determining if a number was prime, but I wonder if either a statistical method or perhaps a sieve would work better. Perhaps it’ll be something to try one day if I get bored.

One additional thing I did was to work out the nested tree structure that you get with something like this. The code is available at GitHub, but here’s a sample run just to whet your appetite:

> (nested-primes #:root 3 #:limit 1000) '(3 (13 (113) (313) (613)) (23 (223) (523) (823)) (43 (443) (643) (743)) (53 (353) (653) (853) (953)) (73 (173) (373) (673) (773)) (83 (283) (383) (683) (883) (983)))

Here’s a nice chart of that:

My goal is to generate a list of all such primes and make a graph of them all, just to see if there’s any sort of visual pattern. I was going to include that in this post, but I have a feeling that code is going to take more than a little time to generate (and actually returning the list directly isn’t probably a good idea as that way it all has to exist in memory at the same time). So I’ll probably optimize this and post it again later. Graph visualization is fun!

If you’d like to download the full code for today’s post, you can do so here: nested primes source

← Confession – Day 51 | All | Confession – Day 52 → |

← Ludum Dare 25 – First favorites | By category | A Brainf**k Interpreter → |

The slowness in your program is caused by the slow prime checker using trial division. A better prime checker, if you don’t mind a statistically tiny chance of error, is a method by Gary Miller and Michael Rabin. See http://programmingpraxis.com/2009/05/01/primality-checking/.

While I can’t help you with a proof, I also suspect the prime checker is most of the issue, especially with bigints. For reference, Paul’s python solution using sympi.isprime takes about 1 second, and my Perl version using Math::Prime::Util::is_prime takes under 0.2s. The former does lots of M-R tests, the latter does a strong BPSW test followed by a BLS75 primality proof.

Thanks both of you.

’tis true that most of the time was spent in the

prime?function. Rerunning it using the solution to the Primality Checking post works much more quickly, although it’s still not quite *that* fast (it took about 15 seconds on the machine I use for testing).That is one nice thing about most batteries-included languages (although not Racket in this case, so far as I know)–you have built in libraries that let you use the nice fast algorithms without having to write them yourself, freeing you to work directly on the problem at hand.

In any case, that code is definitely something to keep in my tool chest. Always nice to have a reasonably fast prime test. And hopefully it will help me generate that nice chart I was talking about earlier. We’ll see.