可截断素数（Truncatable primes）

可截断素数是当从一端连续删除数字时，得到的仍旧是素数的一个素数。

举例

997是左截短素数，因为997, 97 和 7 都是素数。
7393是右截断素数，因为7393, 739, 73 和 7 都是素数。
在可截断质数中不允许有零。

任务

找出小于100万的最大的左可截断和右可截断素数。

相关文章：
寻找给定基数中的最大左可截断素数

C

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
 
#define MAX_PRIME 1000000
char *primes;
int n_primes;
 
/*  Sieve. If we were to handle 10^9 range, use bit field. Regardless,
 *  if a large amount of prime numbers need to be tested, sieve is fast.
 */
void init_primes()
{
	int j;
	primes = malloc(sizeof(char) * MAX_PRIME);
	memset(primes, 1, MAX_PRIME);
	primes[0] = primes[1] = 0;
	int i = 2;
	while (i * i < MAX_PRIME) {
		for (j = i * 2; j < MAX_PRIME; j += i)
			primes[j] = 0;
		while (++i < MAX_PRIME && !primes[i]);
	}
}
 
int left_trunc(int n)
{
	int tens = 1;
	while (tens < n) tens *= 10;
 
	while (n) {
		if (!primes[n]) return 0;
		tens /= 10;
		if (n < tens) return 0;
		n %= tens;
	}
	return 1;
}
 
int right_trunc(int n)
{
	while (n) {
		if (!primes[n]) return 0;
		n /= 10;
	}
	return 1;
}
 
int main()
{
	int n;
	int max_left = 0, max_right = 0;
	init_primes();
 
	for (n = MAX_PRIME - 1; !max_left;  n -= 2)
		if (left_trunc(n)) max_left = n;
 
	for (n = MAX_PRIME - 1; !max_right; n -= 2)
		if (right_trunc(n)) max_right = n;
 
	printf("Left: %d; right: %d\n", max_left, max_right);
	return 0;
}

输出

Left: 998443; right: 739399

对小数(1000000不大)进行质数检验的更快的方法，并自底向上生成可处理的质数：

#include <stdio.h>
 
#define MAXN 1000000
int maxl, maxr;
 
int is_prime(int n)
{
	int p;
	if (n % 3 == 0) return 0;
 
	for (p = 6; p * p <= n; p += 6)
		if (!(n % (p + 1) && n % (p + 5)))
			return 0;
	return 1;
}
 
void left(int n, int tens)
{
	int i, nn;
 
	if (n > maxl) maxl = n;
	if (n < MAXN / 10)
		for (tens *= 10, i = 1; i < 10; i++)
			if (is_prime(nn = i * tens + n))
				left(nn, tens);
}
 
void right(int n)
{
	int i, nn;
	static int d[] = {1,3,7,9};
 
	if (n > maxr) maxr = n;
	if (n < MAXN / 10)
		for (i = 1; i < 4; i++)
			if (is_prime(nn = n * 10 + d[i])) right(nn);
}
 
int main(void)
{
	left(3, 1); left(7, 1);
	right(3); right(5); right(7);
 
	printf("%d %d\n", maxl, maxr);
 
	return 0;
}

输出：

998443 739399

C#

using System;  // 4790@3.6
using System.Collections.Generic;
class truncatable_primes
{
    static void Main()
    {
        uint m = 1000000;
        Console.Write("L " + L(m) + " R " + R(m) + "  ");
        var sw = System.Diagnostics.Stopwatch.StartNew();
        for (int i = 1000; i > 0; i--) { L(m); R(m); }
        Console.Write(sw.Elapsed); Console.Read();
    }
 
    static uint L(uint n)
    {
        n -= n & 1; n--;
        for (uint d, d1 = 100; ; n -= 2)
        {
            while (n % 3 == 0 || n % 5 == 0 || n % 7 == 0) n -= 2;
            if ((d = n % 10) == 3 || d == 7)
            {
                while (d1 < n && d < (d = n % d1) && isP(d)) d1 *= 10;
                if (d1 > n && isP(n)) return n; d1 = 100;
            }
        }
    }
 
    static uint R(uint m)
    {
        var p = new List<uint>() { 2, 3, 5, 7 }; uint n = 20, np;
        for (int i = 1; i < p.Count; n = 10 * p[i++])
        {
            if ((np = n + 1) >= m) break; if (isP(np)) p.Add(np);
            if ((np = n + 3) >= m) break; if (isP(np)) p.Add(np);
            if ((np = n + 7) >= m) break; if (isP(np)) p.Add(np);
            if ((np = n + 9) >= m) break; if (isP(np)) p.Add(np);
        }
        return p[p.Count - 1];
    }
 
    static bool isP(uint n)
    {
        if (n < 7) return n == 2 || n == 3 || n == 5;
        if ((n & 1) == 0 || n % 3 == 0 || n % 5 == 0) return false;
        for (uint r = (uint)Math.Sqrt(n), d = 7; d <= r; d += 30)
            if (n % (d + 00) == 0 || n % (d + 04) == 0 ||
                n % (d + 06) == 0 || n % (d + 10) == 0 ||
                n % (d + 12) == 0 || n % (d + 16) == 0 ||
                n % (d + 22) == 0 || n % (d + 24) == 0) return false;
        return true;
    }
}

输出：

L 998443 R 739399   24 μs

C++

#include <iostream>
#include "prime_sieve.hpp"
 
bool is_left_truncatable(const prime_sieve& sieve, int p) {
    for (int n = 10, q = p; p > n; n *= 10) {
        if (!sieve.is_prime(p % n) || q == p % n)
            return false;
        q = p % n;
    }
    return true;
}
 
bool is_right_truncatable(const prime_sieve& sieve, int p) {
    for (int q = p/10; q > 0; q /= 10) {
        if (!sieve.is_prime(q))
            return false;
    }
    return true;
}
 
int main() {
    const int limit = 1000000;
 
    // find the prime numbers up to the limit
    prime_sieve sieve(limit + 1);
 
    int largest_left = 0;
    int largest_right = 0;
    // find largest left truncatable prime
    for (int p = limit; p >= 2; --p) {
        if (sieve.is_prime(p) && is_left_truncatable(sieve, p)) {
            largest_left = p;
            break;
        }
    }
    // find largest right truncatable prime
    for (int p = limit; p >= 2; --p) {
        if (sieve.is_prime(p) && is_right_truncatable(sieve, p)) {
            largest_right = p;
            break;
        }
    }
    // write results to standard output
    std::cout << "Largest left truncatable prime is " << largest_left << '\n';
    std::cout << "Largest right truncatable prime is " << largest_right << '\n';
    return 0;
}

prime_sieve.hpp

#ifndef PRIME_SIEVE_HPP
#define PRIME_SIEVE_HPP
 
#include <algorithm>
#include <vector>
 
/**
 * A simple implementation of the Sieve of Eratosthenes.
 * See https://en.wikipedia.org/wiki/Sieve_of_Eratosthenes.
 */
class prime_sieve {
public:
    explicit prime_sieve(size_t);
    bool is_prime(size_t) const;
private:
    std::vector<bool> is_prime_;
};
 
/**
 * Constructs a sieve with the given limit.
 *
 * @param limit the maximum integer that can be tested for primality
 */
inline prime_sieve::prime_sieve(size_t limit) {
    limit = std::max(size_t(3), limit);
    is_prime_.resize(limit/2, true);
    for (size_t p = 3; p * p <= limit; p += 2) {
        if (is_prime_[p/2 - 1]) {
            size_t inc = 2 * p;
            for (size_t q = p * p; q <= limit; q += inc)
                is_prime_[q/2 - 1] = false;
        }
    }
}
 
/**
 * Returns true if the given integer is a prime number. The integer
 * must be less than or equal to the limit passed to the constructor.
 *
 * @param n an integer less than or equal to the limit passed to the
 * constructor
 * @return true if the integer is prime
 */
inline bool prime_sieve::is_prime(size_t n) const {
    if (n == 2)
        return true;
    if (n < 2 || n % 2 == 0)
        return false;
    return is_prime_.at(n/2 - 1);
}
 
#endif

输出：

Largest left truncatable prime is 998443
Largest right truncatable prime is 739399

Go

package main
 
import "fmt"
 
func main() {
    sieve(1e6)
    if !search(6, 1e6, "left", func(n, pot int) int { return n % pot }) {
        panic("997?")
    }
    if !search(6, 1e6, "right", func(n, _ int) int { return n / 10 }) {
        panic("7393?")
    }
}
 
var c []bool
 
func sieve(ss int) {
    c = make([]bool, ss)
    c[1] = true
    for p := 2; ; {
        p2 := p * p
        if p2 >= ss {
            break
        }
        for i := p2; i < ss; i += p {
            c[i] = true
        }
        for {
            p++
            if !c[p] {
                break
            }
        }
    }
}
 
func search(digits, pot int, s string, truncFunc func(n, pot int) int) bool {
    n := pot - 1
    pot /= 10
smaller:
    for ; n >= pot; n -= 2 {
        for tn, tp := n, pot; tp > 0; tp /= 10 {
            if tn < tp || c[tn] {
                continue smaller
            }
            tn = truncFunc(tn, tp)
        }
        fmt.Println("max", s, "truncatable:", n)
        return true
    }
    if digits > 1 {
        return search(digits-1, pot, s, truncFunc)
    }
    return false
}

输出：

max left truncatable: 998443
max right truncatable: 739399

Java

import java.util.BitSet;
 
public class Main {
 
	public static void main(String[] args){
 
		final int MAX = 1000000;
 
		//Sieve of Eratosthenes (using BitSet only for odd numbers)
		BitSet primeList = new BitSet(MAX>>1); 
		primeList.set(0,primeList.size(),true); 
 
		int sqroot = (int) Math.sqrt(MAX); 
		primeList.clear(0); 
		for(int num = 3; num <= sqroot; num+=2) 
		{ 
			if( primeList.get(num >> 1) ) 
			{ 
				int inc = num << 1;
				for(int factor = num * num; factor < MAX; factor += inc) 
				{ 
					//if( ((factor) & 1) == 1) 
					//{ 
					primeList.clear(factor >> 1); 
					//} 
				} 
			} 
		}
		//Sieve ends...
 
		//Find Largest Truncatable Prime. (so we start from 1000000 - 1
		int rightTrunc = -1, leftTrunc = -1;
		for(int prime = (MAX - 1) | 1; prime >= 3; prime -= 2)
		{
			if(primeList.get(prime>>1))
			{
				//Already found Right Truncatable Prime?
				if(rightTrunc == -1)
				{
					int right = prime;
					while(right > 0 && right % 2 != 0 && primeList.get(right >> 1)) right /= 10;
					if(right == 0) rightTrunc = prime;
				}
 
				//Already found Left Truncatable Prime?
				if(leftTrunc == -1 )
				{
					//Left Truncation
					String left = Integer.toString(prime);
					if(!left.contains("0"))
					{
						while( left.length() > 0 ){
							int iLeft = Integer.parseInt(left);
							if(!primeList.get( iLeft >> 1)) break;
							left = left.substring(1);
						}
						if(left.length() == 0) leftTrunc = prime;
					}
				}
				if(leftTrunc != -1 && rightTrunc != -1) //Found both? then Stop loop
				{
					break;
				}
			}
		}
		System.out.println("Left  Truncatable : " + leftTrunc);
		System.out.println("Right Truncatable : " + rightTrunc);
	}
}

输出：

Left  Truncatable : 998443
Right Truncatable : 739399

Python

maxprime = 1000000
 
def primes(n):
    multiples = set()
    prime = []
    for i in range(2, n+1):
        if i not in multiples:
            prime.append(i)
            multiples.update(set(range(i*i, n+1, i)))
    return prime
 
def truncatableprime(n):
    'Return a longest left and right truncatable primes below n'
    primelist = [str(x) for x in primes(n)[::-1]]
    primeset = set(primelist)
    for n in primelist:
        # n = 'abc'; [n[i:] for i in range(len(n))] -> ['abc', 'bc', 'c']
        alltruncs = set(n[i:] for i in range(len(n)))
        if alltruncs.issubset(primeset):
            truncateleft = int(n)
            break
    for n in primelist:
        # n = 'abc'; [n[:i+1] for i in range(len(n))] -> ['a', 'ab', 'abc']
        alltruncs = set([n[:i+1] for i in range(len(n))])
        if alltruncs.issubset(primeset):
            truncateright = int(n)
            break
    return truncateleft, truncateright
 
print(truncatableprime(maxprime))

输出：

(998443, 739399)

Ruby

def left_truncatable?(n)
  truncatable?(n) {|i| i.to_s[1..-1].to_i}
end
 
 
def right_truncatable?(n)
  truncatable?(n) {|i| i/10}
end
 
def truncatable?(n, &trunc_func)
  return false if n.to_s.include? "0"
  loop do
    n = trunc_func.call(n)
    return true if n.zero?
    return false unless Prime.prime?(n)
  end
end
 
require 'prime'
primes = Prime.each(1_000_000).to_a.reverse
 
p primes.detect {|p| left_truncatable? p}
p primes.detect {|p| right_truncatable? p}

输出：

998443
739399