Find largest left truncatable prime in a given base: Difference between revisions
(Clarify, I hope!) |
|||
| Line 2: | Line 2: | ||
A [[Truncatable primes|truncatable prime]] is one where all non-empty substrings that finish at the end of the number (right-substrings) are also primes ''when understood as numbers in a particular base''. The largest such prime in a given (integer) base is therefore computable, provided the base is larger than 2. |
A [[Truncatable primes|truncatable prime]] is one where all non-empty substrings that finish at the end of the number (right-substrings) are also primes ''when understood as numbers in a particular base''. The largest such prime in a given (integer) base is therefore computable, provided the base is larger than 2. |
||
Let's consider what happens in base 10. Obviously the right most digit must be prime, so in base 10 candidates are 2,3,5,7. Putting a digit in the range 1 to base-1 in front of each candidate must result in a prime. So 2 and 5, like the whale and the petunias in ''The Hitchhiker's Guide to the Galaxy'', come into |
Let's consider what happens in base 10. Obviously the right most digit must be prime, so in base 10 candidates are 2,3,5,7. Putting a digit in the range 1 to base-1 in front of each candidate must result in a prime. So 2 and 5, like the whale and the petunias in ''The Hitchhiker's Guide to the Galaxy'', come into existence only to be extinguished before they have time to realize it, because 2 and 5 preceded by any digit in the range 1 to base-1 is not prime. Some numbers formed by preceding 3 or 7 by a digit in the range 1 to base-1 are prime. So 13,17,23,37,43,47,53,67,73,83,97 are candidates. Again, putting a digit in the range 1 to base-1 in front of each candidate must be a prime. Repeating until there are no larger candidates finds the largest left truncatable prime. |
||
Let's work base 3 by hand: |
Let's work base 3 by hand: |
||
Revision as of 12:29, 3 September 2013
You are encouraged to solve this task according to the task description, using any language you may know.
A truncatable prime is one where all non-empty substrings that finish at the end of the number (right-substrings) are also primes when understood as numbers in a particular base. The largest such prime in a given (integer) base is therefore computable, provided the base is larger than 2.
Let's consider what happens in base 10. Obviously the right most digit must be prime, so in base 10 candidates are 2,3,5,7. Putting a digit in the range 1 to base-1 in front of each candidate must result in a prime. So 2 and 5, like the whale and the petunias in The Hitchhiker's Guide to the Galaxy, come into existence only to be extinguished before they have time to realize it, because 2 and 5 preceded by any digit in the range 1 to base-1 is not prime. Some numbers formed by preceding 3 or 7 by a digit in the range 1 to base-1 are prime. So 13,17,23,37,43,47,53,67,73,83,97 are candidates. Again, putting a digit in the range 1 to base-1 in front of each candidate must be a prime. Repeating until there are no larger candidates finds the largest left truncatable prime.
Let's work base 3 by hand:
0 and 1 are not prime so the last digit must be 2. 123 = 510 which is prime, 223 = 810 which is not so 123 is the only candidate. 1123 = 1410 which is not prime, 2123 = 2310 which is, so 2123 is the only candidate. 12123 = 5010 which is not prime, 22123 = 7710 which also is not prime. So there are no more candidates, therefore 23 is the largest left truncatable prime in base 3.
The task is to reconstruct as much, and possibly more, of the table in the OEIS as you are able.
Related Tasks:
BBC BASIC
Uses the Huge Integer Math & Encryption library from http://devotechs.com/ <lang bbcbasic> HIMEM = PAGE + 3000000
INSTALL @lib$+"HIMELIB"
PROC_himeinit("HIMEkey")
DIM old$(20000), new$(20000)
h1% = 1 : h2% = 2 : h3% = 3 : h4% = 4
FOR base% = 3 TO 17
PRINT "Base "; base% " : " FN_largest_left_truncated_prime(base%)
NEXT
END
DEF FN_largest_left_truncated_prime(base%)
LOCAL digit%, i%, new%, old%, prime%, fast%, slow%
fast% = 1 : slow% = 50
old$() = ""
PROC_hiputdec(1, STR$(base%))
PROC_hiputdec(2, "1")
REPEAT
new% = 0 : new$() = ""
PROC_hiputdec(3, "0")
FOR digit% = 1 TO base%-1
SYS `hi_Add`, ^h2%, ^h3%, ^h3%
FOR i% = 0 TO old%-1
PROC_hiputdec(4, old$(i%))
SYS `hi_Add`, ^h3%, ^h4%, ^h4%
IF old% OR digit% > 1 THEN
IF old% > 100 THEN
SYS `hi_IsPrime_RB`, ^fast%, ^h4% TO prime%
ELSE
SYS `hi_IsPrime_RB`, ^slow%, ^h4% TO prime%
ENDIF
IF prime% THEN new$(new%) = FN_higetdec(4) : new% += 1
ENDIF
NEXT
NEXT
SYS `hi_Mul`, ^h1%, ^h2%, ^h2%
SWAP old$(), new$()
SWAP old%, new%
UNTIL old% = 0
= new$(new%-1)</lang>
Output:
Base 3 : 23 Base 4 : 4091 Base 5 : 7817 Base 6 : 4836525320399 Base 7 : 817337 Base 8 : 14005650767869 Base 9 : 1676456897 Base 10 : 357686312646216567629137 Base 11 : 2276005673 Base 12 : 13092430647736190817303130065827539 Base 13 : 812751503 Base 14 : 615419590422100474355767356763 Base 15 : 34068645705927662447286191 Base 16 : 1088303707153521644968345559987 Base 17 : 13563641583101
C
<lang c>#include <stdio.h>
- include <gmp.h>
typedef unsigned long ulong;
ulong small_primes[] = {2,3,5,7,11,13,17,19,23,29,31,37,41, 43,47,53,59,61,67,71,73,79,83,89,97};
- define MAX_STACK 128
mpz_t tens[MAX_STACK], value[MAX_STACK], answer;
ulong base, seen_depth;
void add_digit(ulong i) { ulong d; for (d = 1; d < base; d++) { mpz_set(value[i], value[i-1]); mpz_addmul_ui(value[i], tens[i], d); if (!mpz_probab_prime_p(value[i], 1)) continue;
if (i > seen_depth || (i == seen_depth && mpz_cmp(value[i], answer) == 1)) { if (!mpz_probab_prime_p(value[i], 50)) continue;
mpz_set(answer, value[i]); seen_depth = i; gmp_fprintf(stderr, "\tb=%lu d=%2lu | %Zd\n", base, i, answer); }
add_digit(i+1); } }
void do_base() { ulong i; mpz_set_ui(answer, 0); mpz_set_ui(tens[0], 1); for (i = 1; i < MAX_STACK; i++) mpz_mul_ui(tens[i], tens[i-1], base);
for (seen_depth = i = 0; small_primes[i] < base; i++) { fprintf(stderr, "\tb=%lu digit %lu\n", base, small_primes[i]); mpz_set_ui(value[0], small_primes[i]); add_digit(1); } gmp_printf("%d: %Zd\n", base, answer); }
int main(void) { ulong i; for (i = 0; i < MAX_STACK; i++) { mpz_init_set_ui(tens[i], 0); mpz_init_set_ui(value[i], 0); } mpz_init_set_ui(answer, 0);
for (base = 22; base < 30; base++) do_base();
return 0; }</lang>
- Output:
$ ./a.out 2&>/dev/null 22: 33389741556593821170176571348673618833349516314271 23: 116516557991412919458949 24: 10594160686143126162708955915379656211582267119948391137176997290182218433 25: 8211352191239976819943978913 26: 12399758424125504545829668298375903748028704243943878467 27: 10681632250257028944950166363832301357693 28: 720639908748666454129676051084863753107043032053999738835994276213 29: 4300289072819254369986567661 ...
Haskell
Miller-Rabin test code from HaskellWiki, with modifications. <lang haskell>primesTo100 = [2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61,67,71,73,79,83,89,97]
-- (eq. to) find2km (2^k * n) = (k,n) find2km :: Integral a => a -> (Int,a) find2km n = f 0 n where f k m
| r == 1 = (k,m)
| otherwise = f (k+1) q
where (q,r) = quotRem m 2
-- n is the number to test; a is the (presumably randomly chosen) witness millerRabinPrimality :: Integer -> Integer -> Bool millerRabinPrimality n a
| a >= n_ = True
| b0 == 1 || b0 == n_ = True
| otherwise = iter (tail b)
where
n_ = n-1
(k,m) = find2km n_
b0 = powMod n a m
b = take k $ iterate (squareMod n) b0
iter [] = False
iter (x:xs)
| x == 1 = False
| x == n_ = True
| otherwise = iter xs
-- (eq. to) pow_ (*) (^2) n k = n^k pow_ :: (Num a, Integral b) => (a->a->a) -> (a->a) -> a -> b -> a pow_ _ _ _ 0 = 1 pow_ mul sq x_ n_ = f x_ n_ 1
where
f x n y
| n == 1 = x `mul` y
| r == 0 = f x2 q y
| otherwise = f x2 q (x `mul` y)
where
(q,r) = quotRem n 2
x2 = sq x
mulMod :: Integral a => a -> a -> a -> a mulMod a b c = (b * c) `mod` a squareMod :: Integral a => a -> a -> a squareMod a b = (b * b) `rem` a
-- (eq. to) powMod m n k = n^k `mod` m powMod :: Integral a => a -> a -> a -> a powMod m = pow_ (mulMod m) (squareMod m)
-- Caller supplies a witness list w, which may be used for MR test. -- Use faster trial division against a small primes list first, to -- weed out more obvious composites. is_prime w n | n < 100 = n `elem` primesTo100 | any ((==0).(n`mod`)) primesTo100 = False | otherwise = all (millerRabinPrimality n) w
-- final result gets a more thorough Miller-Rabin left_trunc base = head $ filter (is_prime primesTo100) (reverse hopeful) where hopeful = extend base $ takeWhile (<base) primesTo100 where extend b x = if null d then x else extend (b*base) d where d = concatMap addDigit [1..base-1] -- we do *one* prime test, which seems good enough in practice addDigit a = filter (is_prime [3]) $ map (a*b+) x
main = mapM_ print $ map (\x->(x, left_trunc x)) [3..21]</lang>
(3,23) (4,4091) (5,7817) (6,4836525320399) (7,817337) (8,14005650767869) (9,1676456897) (10,357686312646216567629137) (11,2276005673) (12,13092430647736190817303130065827539) (13,812751503) (14,615419590422100474355767356763) (15,34068645705927662447286191) (16,1088303707153521644968345559987) (17,13563641583101) (18,571933398724668544269594979167602382822769202133808087) (19,546207129080421139) (20,1073289911449776273800623217566610940096241078373) (21,391461911766647707547123429659688417)
J
I'm fairly sure this is correct but it sure has long run time for base 12. <lang J>Filter=: (#~`)(`:6) Filter=: (#~`)(`:6)
f=: 3 : 0 BASE=. y primes=. (1&p:)Filter DIGITS=. }. i. x: BASE R=. 0 A=. ,0 while. #A do.
R=. >: R
B=. A
A=. primes , , A +/ BASE #. (_,R) {. ,. DIGITS
end. >./ B ) l NB. f N where N is the base.
(,. f"0) (3+i.8) 3 23 4 4091 5 7817 6 4836525320399 7 817337 8 14005650767869 9 1676456897
10 357686312646216567629137
f 11
2276005673
</lang> What can I say? J is an executable mathematical notation better than algebraic notation. That is the point, not that the execution time of a particular implementation is such and such. The current implementation of J uses custom arbitrary length integer code instead of gmp (j is not optimized to the extent gmp is optimized with its assembly code for various platforms. J is now open source. Join the project and install gmp!). J is terser than algebra. J handles text. J is simpler than algebra. J's precedence rules of noun-verb sequences evaluate strictly right to left. Versus the 17 precedence rules of C++. Versus algebraic confusion. "What is the intent of "force/area*time"? Surprise! It probably means physical power, force/(area*time). Different precedence, confusing, from the algebraic evaluation of, in c or mathematica, 8/3*7 == 7*8/3. Consider the simplifications of electrodynamics equations from Maxwell's original 20 to the 4 in modern notation. J benefits thought. Think in J. www.jsoftware.com Clearly, the various ways to express a thought in J perform differently upon execution. Languages like forth, python, c were developed AS PROGRAMMING LANGUAGES. J was developed as a TOOL OF THOUGHT. To untrained eye J looks like line noise with bit 8 clear. Get over it. Colon (:) and dot (.) are inflections, creating related concepts to the uninflected idea. Dyadic + means add. Inflected dyadic +, means "greatest common denominator" which for binary operands is "or". Am I wrong? Does the j solution not have 3 times fewer lines of code than the other languages have? Are my programs typically longer than those programs of others? Honeycombs/Python Honeycombs#C
Java
Code:
<lang java>import java.math.BigInteger; import java.util.*;
class LeftTruncatablePrime {
private static List<BigInteger> getNextLeftTruncatablePrimes(BigInteger n, int radix, int millerRabinCertainty)
{
List<BigInteger> probablePrimes = new ArrayList<BigInteger>();
String baseString = n.equals(BigInteger.ZERO) ? "" : n.toString(radix);
for (int i = 1; i < radix; i++)
{
BigInteger p = new BigInteger(Integer.toString(i, radix) + baseString, radix);
if (p.isProbablePrime(millerRabinCertainty))
probablePrimes.add(p);
}
return probablePrimes;
}
public static BigInteger getLargestLeftTruncatablePrime(int radix, int millerRabinCertainty)
{
List<BigInteger> lastList = null;
List<BigInteger> list = getNextLeftTruncatablePrimes(BigInteger.ZERO, radix, millerRabinCertainty);
while (!list.isEmpty())
{
lastList = list;
list = new ArrayList<BigInteger>();
for (BigInteger n : lastList)
list.addAll(getNextLeftTruncatablePrimes(n, radix, millerRabinCertainty));
}
if (lastList == null)
return null;
Collections.sort(lastList);
return lastList.get(lastList.size() - 1);
}
public static void main(String[] args)
{
int maxRadix = Integer.parseInt(args[0]);
int millerRabinCertainty = Integer.parseInt(args[1]);
for (int radix = 3; radix <= maxRadix; radix++)
{
BigInteger largest = getLargestLeftTruncatablePrime(radix, millerRabinCertainty);
System.out.print("n=" + radix + ": ");
if (largest == null)
System.out.println("No left-truncatable prime");
else
System.out.println(largest + " (in base " + radix + "): " + largest.toString(radix));
}
}
}</lang>
Example:
java LeftTruncatablePrime 17 1000000 n=3: 23 (in base 3): 212 n=4: 4091 (in base 4): 333323 n=5: 7817 (in base 5): 222232 n=6: 4836525320399 (in base 6): 14141511414451435 n=7: 817337 (in base 7): 6642623 n=8: 14005650767869 (in base 8): 313636165537775 n=9: 1676456897 (in base 9): 4284484465 n=10: 357686312646216567629137 (in base 10): 357686312646216567629137 n=11: 2276005673 (in base 11): a68822827 n=12: 13092430647736190817303130065827539 (in base 12): 471a34a164259ba16b324ab8a32b7817 n=13: 812751503 (in base 13): cc4c8c65 n=14: 615419590422100474355767356763 (in base 14): d967ccd63388522619883a7d23 n=15: 34068645705927662447286191 (in base 15): 6c6c2ce2ceeea4826e642b n=16: 1088303707153521644968345559987 (in base 16): dbc7fba24fe6aec462abf63b3 n=17: 13563641583101 (in base 17): 6c66cc4cc83
Mathematica
<lang>LargestLeftTruncatablePrimeInBase[n_] :=
Max[NestWhile[{Select[
Flatten@Outer[Function[{a, b}, #2 a + b],
Range[1, n - 1], #1], PrimeQ], n #2} &, {{0},
1}, #1 != {} &, 1, Infinity, -1]1]</lang>
Example:
<lang>Do[Print[n, "\t", LargestLeftTruncatablePrimeInBase@n], {n, 3, 17}]</lang>
Output:
3 23 4 4091 5 7817 6 4836525320399 7 817337 8 14005650767869 9 1676456897 10 357686312646216567629137 11 2276005673 12 13092430647736190817303130065827539 13 812751503 14 615419590422100474355767356763 15 34068645705927662447286191 16 1088303707153521644968345559987 17 13563641583101
PARI/GP
Takes half a second to find the terms up to 10, with proofs of primality. The time can be halved without proofs (use ispseudoprime in place of isprime).
<lang parigp>a(n)=my(v=primes(primepi(n-1)),u,t,b=1,best); while(#v, best=vecmax(v); b*=n; u=List(); for(i=1,#v,for(k=1,n-1,if(isprime(t=v[i]+k*b), listput(u,t)))); v=Vec(u)); best</lang>
Perl 6
Using the Miller-Rabin definition of is-prime from Miller-Rabin primality test#Perl 6, or the built-in, if available. For speed we do a single try, which allows a smattering of composites in, but this does not appear to damage the algorithm much, and the correct answer appears to be within a few candidates of the end all the time. We keep the last few hundred candidates, sort them into descending order, then pick the largest that passes Miller-Rabin with 100 tries.
I believe the composites don't hurt the algorithm because they do not confer any selective advantage on their "children", so their average length is no longer than the correct solution. In addition, the candidates in the finalist set all appear to be largely independent of each other, so any given composite tends to contribute only a single candidate, if any. I have no proof of this; I prefer my math to have more vigor than rigor. :-) <lang perl6>for 3..* -> $base {
say "Starting base $base...";
my @stems = grep { is-prime($_, 100)}, ^$base;
my @runoff;
for 1 .. * -> $digits {
print ' ', @stems.elems;
my @new;
my $place = $base ** $digits;
my $tries = 1;
for @stems -> $stem {
for 1 ..^ $base -> $digit { my $new = $digit * $place + $stem; @new.push($new) if is-prime($new, $tries);
}
}
last unless @new;
push @runoff, @stems if @new < @stems and @new < 100;
@stems = @new;
}
push @runoff, @stems;
say "\n Finalists: ", +@runoff;
for @runoff.sort(-*) -> $finalist {
my $b = $finalist.base($base); say " Checking :$base\<", $b, '>'; my $ok = True; for $base ** 2, $base ** 3, $base ** 4 ... $base ** $b.chars -> $place { my $f = $finalist % $place; if not is-prime($f, 100) { say " Oops, :$base\<", $f.base($base), '> is not a prime!!!'; $ok = False; last; } } next unless $ok;
say " Largest ltp in base $base = $finalist"; last;
}
}</lang>
- Output:
Starting base 3...
1 1 1
Finalists: 1
Checking :3<212>
Largest ltp in base 3 = 23
Starting base 4...
2 3 5 5 6 4 2
Finalists: 12
Checking :4<3323233>
Oops, :4<33> is not a prime!!!
Checking :4<1323233>
Oops, :4<33> is not a prime!!!
Checking :4<333323>
Largest ltp in base 4 = 4091
Starting base 5...
2 4 4 3 1 1
Finalists: 8
Checking :5<222232>
Largest ltp in base 5 = 7817
Starting base 6...
3 4 12 26 45 56 65 67 66 57 40 25 14 9 4 2 1
Finalists: 285
Checking :6<14141511414451435>
Largest ltp in base 6 = 4836525320399
Starting base 7...
3 6 6 4 1 1 1
Finalists: 11
Checking :7<6642623>
Largest ltp in base 7 = 817337
Starting base 8...
4 12 29 50 67 79 65 52 39 28 17 8 3 2 1
Finalists: 294
Checking :8<313636165537775>
Largest ltp in base 8 = 14005650767869
Starting base 9...
4 9 15 17 24 16 9 6 5 3
Finalists: 63
Checking :9<4284484465>
Largest ltp in base 9 = 1676456897
Starting base 10...
4 11 40 101 197 335 439 536 556 528 456 358 278 212 117 72 42 24 13 6 5 4 3 1
Finalists: 287
Checking :10<357686312646216567629137>
Largest ltp in base 10 = 357686312646216567629137
Starting base 11...
4 8 15 18 15 8 4 2 1
Finalists: 48
Checking :11<A68822827>
Largest ltp in base 11 = 2276005673
Starting base 12...
5 24 124 431 1179 2616 4948 8054 11732 15460 18100 19546 19777 18280 15771 12574 9636 6866 4625 2990 1867 1152 627 357 189 107 50 20 9 3 1 1
Finalists: 190
Checking :12<471A34A164259BA16B324AB8A32B7817>
Largest ltp in base 12 = 13092430647736190817303130065827539
Starting base 13...
5 13 21 23 17 11 7 4
Finalists: 62
Checking :13<CC4C8C65>
Largest ltp in base 13 = 812751503
Starting base 14...
6 28 103 308 694 1329 2148 3089 3844 4290 4311 3923 3290 2609 1948 1298 809 529 332 167 97 55 27 13 5 2
Finalists: 366
Checking :14<D967CCD63388522619883A7D23>
Largest ltp in base 14 = 615419590422100474355767356763
Starting base 15...
6 22 80 213 401 658 955 1208 1307 1209 1075 841 626 434 268 144 75 46 27 14 7 1
Finalists: 314
Checking :15<6C6C2CE2CEEEA4826E642B>
Largest ltp in base 15 = 34068645705927662447286191
Starting base 16...
6 32 131 355 792 1369 2093 2862 3405 3693 3519 3140 2660 1930 1342 910 571 335 180 95 50 25 9 5 1
Finalists: 365
Checking :16<DBC7FBA24FE6AEC462ABF63B3>
Largest ltp in base 16 = 1088303707153521644968345559987
Starting base 17...
6 23 47 64 80 62 43 32 23 8 1
Finalists: 249
Checking :17<6C66CC4CC83>
Largest ltp in base 17 = 13563641583101
Starting base 18...
7 49 311 1396 5117 15243 38818 85684 167133 293518 468456 680171 911723 1133959 1313343 1423597 1449405 1395514 1270222 1097353 902477 707896 529887 381239 263275 174684 111046 67969 40704 23201 12793 6722 3444 1714 859 422 205 98 39 14 7 1 1
Finalists: 364
Checking :18<AF93E41A586HE75A7HHAAB7HE12FG79992GA7741B3D>
Largest ltp in base 18 = 571933398724668544269594979167602382822769202133808087
Starting base 19...
7 29 61 106 122 117 93 66 36 18 10 10 6 4
Finalists: 350
Checking :19<CIEG86GCEA2C6H>
Largest ltp in base 19 = 546207129080421139
Starting base 20...
8 56 321 1311 4156 10963 24589 47737 83011 129098 181707 234685 278792 306852 315113 302684 273080 233070 188331 145016 105557 73276 48819 31244 19237 11209 6209 3383 1689 917 430 196 80 44 20 7 2
Finalists: 349
Checking :20<FC777G3CG1FIDI9I31IE5FFB379C7A3F6EFID>
Largest ltp in base 20 = 1073289911449776273800623217566610940096241078373
Starting base 21...
8 41 165 457 1079 2072 3316 4727 6003 6801 7051 6732 5862 4721 3505 2474 1662 1039 628 369 219 112 52 17 13 4 2
Finalists: 200
Checking :21<G8AGG2GCA8CAK4K68GEA4G2K22H>
Largest ltp in base 21 = 391461911766647707547123429659688417
Starting base 22...
8 48 261 1035 3259 8247 17727 33303 55355 82380 111919 137859 158048 167447 165789 153003 132516 108086 83974 61950 43453 29212 18493 11352 6693 3738 2053 1125 594 293 145 70 31 13 6 2 1
Finalists: 268
Checking :22<FFHALC8JFB9JKA2AH9FAB4I9L5I9L3GF8D5L5>
Largest ltp in base 22 = 33389741556593821170176571348673618833349516314271
Starting base 23...
8 30 78 137 181 200 186 171 121 100 67 41 24 16 9 2 1
Finalists: 260
Checking :23<IMMGM6C6IMCI66A4H>
Largest ltp in base 23 = 116516557991412919458949
Python
<lang python>import random
def is_probable_prime(n,k):
#this uses the miller-rabin primality test found from rosetta code
if n==0 or n==1:
return False
if n==2:
return True
if n % 2 == 0:
return False
s = 0
d = n-1
while True:
quotient, remainder = divmod(d, 2)
if remainder == 1:
break
s += 1
d = quotient
def try_composite(a):
if pow(a, d, n) == 1:
return False
for i in range(s):
if pow(a, 2**i * d, n) == n-1:
return False
return True # n is definitely composite
for i in range(k):
a = random.randrange(2, n)
if try_composite(a):
return False
return True # no base tested showed n as composite
def largest_left_truncatable_prime(base):
radix = 0
candidates = [0]
while True:
new_candidates=[]
multiplier = base**radix
for i in range(1,base):
new_candidates += [x+i*multiplier for x in candidates if is_probable_prime(x+i*multiplier,30)]
if len(new_candidates)==0:
return max(candidates)
candidates = new_candidates
radix += 1
for b in range(3,24):
print("%d:%d\n" % (b,largest_left_truncatable_prime(b)))
</lang>
Output:
3:23 4:4091 5:7817 6:4836525320399 7:817337 8:14005650767869 9:1676456897 10:357686312646216567629137 11:2276005673 12:13092430647736190817303130065827539 13:812751503
Ruby
Ruby Ruby
<lang ruby>
- Compute the largest left truncatable prime
- Nigel_Galloway
- September 15th., 2012.
require 'prime' BASE = 3 MAX = 500 stems = Prime.each(BASE-1).to_a (1..MAX-1).each {|i|
print "#{stems.length} "
t = []
b = BASE ** i
stems.each {|z|
(1..BASE-1).each {|n|
c = n*b+z
t.push(c) if c.prime?
}}
break if t.empty?
stems = t
} puts "The largest left truncatable prime #{"less than #{BASE ** MAX} " if MAX < 500}in base #{BASE} is #{stems.max}" </lang> By changing BASE from 3 to 14 this produces the solutions in 'Number of left truncatable primes in a given base' on the Discussion Page for bases except 10, 12 and 14.
The maximum left truncatable prime in bases 10 , 12, and 14 are very large. By changing MAX to 6 and BASE to 10 solves related task 1:
The largest left truncatable prime less than 1000000 in base 10 is 998443
JRuby
I require a fast probably prime test. Java has one, is it any good? Let's find out. Ruby can borrow from Java using JRuby. Modifying the Ruby solution: <lang Ruby>
- Compute the largest left truncatable prime
- Nigel_Galloway
- September 15th., 2012.
require 'prime' require 'java' BASE = 10 MAX = 500 stems = Prime.each(BASE-1).to_a (1..MAX-1).each {|i|
print "#{stems.length} "
t = []
b = BASE ** i
stems.each {|z|
(1..BASE-1).each {|n|
c = n*b+z
t.push(c) if java.math.BigInteger.new(c.to_s).isProbablePrime(100)
}}
break if t.empty?
stems = t
} puts "\nThe largest left truncatable prime #{"less than #{BASE ** MAX} " if MAX < 500}in base #{BASE} is #{stems.max}" </lang> Produces all the reults in 'Number of left truncatable primes in a given base' on the discussion page. For bases 18, 20, and 22 I changed the confidence level from 100 to 5 and checked the final answer. Even so base 18 takes a while. For base 24:
9 87 677 3808 17096 63509 199432 545332 1319708 2863180 Error: Your application used more memory than the safety cap of 500m. Specify -J-Xmx####m to increase it (#### = cap size in MB). Specify -w for full OutOfMemoryError stack trace
That is going to be big!
Racket
<lang racket>
- lang racket
(require math/number-theory)
(define (prepend-digit b d i n)
(+ (* d (expt b i)) n))
(define (extend b i ts)
(define ts*
(for/list ([t (in-set ts)])
(for/set ([d (in-range 1 b)]
#:when (prime? (prepend-digit b d i t)))
(prepend-digit b d i t))))
(apply set-union ts*))
(define (truncables b n)
; return set of truncables of length n in base b
(if (= n 1)
(for/set ([d (in-range 1 b)] #:when (prime? d)) d)
(extend b (- n 1) (truncables b (- n 1)))))
(define (largest b)
(let loop ([ts (truncables b 1)]
[n 1])
(define ts* (extend b n ts))
(if (set-empty? ts*)
(apply max (set->list ts))
(loop ts* (+ n 1)))))
(for/list ([b (in-range 3 18)])
(define l (largest b)) ; (displayln (list b l)) (list b l))
- Output
'((3 23)
(4 4091) (5 7817) (6 4836525320399) (7 817337) (8 14005650767869) (9 1676456897) (10 357686312646216567629137) (11 2276005673) (12 13092430647736190817303130065827539) (13 812751503) (14 615419590422100474355767356763) (15 34068645705927662447286191) (16 1088303707153521644968345559987) (17 13563641583101))
</lang>
Tcl
<lang tcl>package require Tcl 8.5
proc tcl::mathfunc::modexp {a b n} {
for {set c 1} {$b} {set a [expr {$a*$a%$n}]} {
if {$b & 1} {
set c [expr {$c*$a%$n}]
}
set b [expr {$b >> 1}]
}
return $c
}
- Based on Miller-Rabin primality testing, but with small prime check first
proc is_prime {n {count 10}} {
# fast check against small primes
foreach p {
2 3 5 7 11 13 17 19 23 29 31 37 41 43 47 53 59 61 67 71 73 79 83 89 97
} {
if {$n == $p} {return true} if {$n % $p == 0} {return false}
}
# write n-1 as 2^s·d with d odd by factoring powers of 2 from n-1
set d [expr {$n - 1}]
for {set s 0} {$d & 1 == 0} {incr s} {
set d [expr {$d >> 1}]
}
for {} {$count > 0} {incr count -1} {
set a [expr {2 + int(rand()*($n - 4))}]
set x [expr {modexp($a, $d, $n)}]
if {$x == 1 || $x == $n - 1} continue
for {set r 1} {$r < $s} {incr r} {
set x [expr {modexp($x, 2, $n)}]
if {$x == 1} {return false}
if {$x == $n - 1} break
}
if {$x != $n-1} {return false}
} return true
}
proc max_left_truncatable_prime {base} {
set stems {}
for {set i 2} {$i < $base} {incr i} {
if {[is_prime $i]} { lappend stems $i }
}
set primes $stems
set size 0
for {set b $base} {[llength $stems]} {set b [expr {$b * $base}]} {
# Progress monitoring is nice once we get to 10 and beyond... if {$base > 9} { puts "\t[llength $stems] candidates at length [incr size]" } set primes $stems set certainty [expr {[llength $primes] > 100 ? 1 : 5}] set stems {} foreach s $primes { for {set i 1} {$i < $base} {incr i} { set n [expr {$b*$i + $s}] if {[is_prime $n $certainty]} { lappend stems $n } } }
}
# Could be several at same length; choose largest
return [tcl::mathfunc::max {*}$primes]
}
for {set i 3} {$i <= 20} {incr i} {
puts "$i: [max_left_truncatable_prime $i]"
}</lang>
- Output up to base 12 (tab-indented parts are progress messages):
3: 23 4: 4091 5: 7817 6: 4836525320399 7: 817337 8: 14005650767869 9: 1676456897 4 candidates at length 1 11 candidates at length 2 39 candidates at length 3 99 candidates at length 4 192 candidates at length 5 326 candidates at length 6 429 candidates at length 7 521 candidates at length 8 545 candidates at length 9 517 candidates at length 10 448 candidates at length 11 354 candidates at length 12 276 candidates at length 13 212 candidates at length 14 117 candidates at length 15 72 candidates at length 16 42 candidates at length 17 24 candidates at length 18 13 candidates at length 19 6 candidates at length 20 5 candidates at length 21 4 candidates at length 22 3 candidates at length 23 1 candidates at length 24 10: 357686312646216567629137 4 candidates at length 1 8 candidates at length 2 15 candidates at length 3 18 candidates at length 4 15 candidates at length 5 8 candidates at length 6 4 candidates at length 7 2 candidates at length 8 1 candidates at length 9 11: 2276005673 5 candidates at length 1 23 candidates at length 2 119 candidates at length 3 409 candidates at length 4 1126 candidates at length 5 2504 candidates at length 6 4746 candidates at length 7 7727 candidates at length 8 11257 candidates at length 9 14860 candidates at length 10 17375 candidates at length 11 18817 candidates at length 12 19027 candidates at length 13 17594 candidates at length 14 15192 candidates at length 15 12106 candidates at length 16 9292 candidates at length 17 6621 candidates at length 18 4466 candidates at length 19 2889 candidates at length 20 1799 candidates at length 21 1109 candidates at length 22 601 candidates at length 23 346 candidates at length 24 181 candidates at length 25 103 candidates at length 26 49 candidates at length 27 19 candidates at length 28 8 candidates at length 29 2 candidates at length 30 1 candidates at length 31 1 candidates at length 32 12: 13092430647736190817303130065827539 5 candidates at length 1 13 candidates at length 2 20 candidates at length 3 23 candidates at length 4 17 candidates at length 5 11 candidates at length 6 7 candidates at length 7 4 candidates at length 8 13: 812751503
I think I'll need to find a faster computer to calculate much more of the sequence, but memory consumption is currently negligible so there's no reason to expect there to be any major problems.