Carmichael 3 strong pseudoprimes
You are encouraged to solve this task according to the task description, using any language you may know.
A lot of composite numbers can be detected by the Miller Rabin Test, but there are some that evade it.
The purpose of this task is to investigate such numbers. The method suggested is based on Notes by G.J.O Jameson March 2010
The objective is to find Carmichael numbers of the form Prime1 X Prime2 X Prime3.
- Prime1 < Prime2 < Prime3 for all Prime1 upto 61 see page 7 of Notes by G.J.O Jameson March 2010.
For a given Prime1
for 1 < h3 < Prime1
- g=h3+Prime1
- for 0 < d < h3+Prime1
- if (h3+Prime1)*(Prime1-1) mod d == 0 and -Prime1 squared mod h3 == d mod h3
- then
- Prime2 = 1 + ((Prime1-1) * g/d)
- next d if Prime2 is not prime
- Prime3 = 1 + (Prime1*Prime2/h3)
- next d if Prime3 is not prime
- next d if (Prime2*Prime3) mod (Prime1-1) not equal 1
- Prime1 * Prime2 * Prime3 is a Carmichael Number
Related Tasks:
D
From the Python entry, with several changes. <lang d>import std.stdio, std.algorithm, std.range, std.bitmanip;
/// Extensible sieve of Eratosthenes. struct IsPrime {
//__gshared static BitArray multiples = [true, true, false]; __gshared static BitArray multiples;
static this() {
multiples.length = 3;
multiples[0] = true;
multiples[1] = true;
}
static bool opCall(in size_t n) /*nothrow*/ {
if (n >= multiples.length) {
// Extend the sieve.
immutable oldLen = multiples.length; // Not nothrow.
immutable newMax = max((oldLen - 1) * 2, n);
multiples.length = newMax + 1;
foreach (immutable p; 2 .. oldLen)
if (!multiples[p])
for (size_t i = p * ((oldLen + p) / p);
i < newMax + 1; i += p)
multiples[i] = true;
foreach (immutable i; oldLen .. newMax + 1)
if (!multiples[i])
for (size_t j = i * 2; j < newMax + 1; j += i)
multiples[j] = true;
}
return !multiples[n];
}
}
struct Carmichael {
immutable int p1;
static int mod(in int n, in int m) pure nothrow {
return ((n % m) + m) % m;
}
int opApply(immutable int delegate(in ref int[3]) dg) {
int result;
if (IsPrime(p1))
foreach (immutable h3; 2 .. p1) {
immutable g = h3 + p1;
foreach (immutable d; 1 .. g)
if ((g * (p1 - 1)) % d == 0
&& mod(-p1 * p1, h3) == d % h3) {
immutable p2 = 1 + ((p1 - 1) * g / d);
if (IsPrime(p2)) {
immutable p3 = 1 + (p1 * p2 / h3);
if (IsPrime(p3))
if ((p2 * p3) % (p1 - 1) == 1) {
// Heap-allocates:
//int[3] triple = [p1, p2, p3];
int[3] triple = void;
triple[0] = p1;
triple[1] = p2;
triple[2] = p3;
result = dg(triple);
if (result) break;
}
}
}
}
return result; }
}
void main() {
foreach (immutable n; 0 .. 62)
if (IsPrime(n))
foreach (const c; Carmichael(n))
writefln("%(%d x %)", c);
}</lang>
- Output:
3 x 11 x 17 5 x 29 x 73 5 x 17 x 29 5 x 13 x 17 7 x 19 x 67 7 x 31 x 73 7 x 13 x 31 7 x 23 x 41 7 x 73 x 103 7 x 13 x 19 13 x 61 x 397 13 x 37 x 241 13 x 97 x 421 13 x 37 x 97 13 x 37 x 61 17 x 41 x 233 17 x 353 x 1201 19 x 43 x 409 19 x 199 x 271 23 x 199 x 353 29 x 113 x 1093 29 x 197 x 953 31 x 991 x 15361 31 x 61 x 631 31 x 151 x 1171 31 x 61 x 271 31 x 61 x 211 31 x 271 x 601 31 x 181 x 331 37 x 109 x 2017 37 x 73 x 541 37 x 613 x 1621 37 x 73 x 181 37 x 73 x 109 41 x 1721 x 35281 41 x 881 x 12041 41 x 101 x 461 41 x 241 x 761 41 x 241 x 521 41 x 73 x 137 41 x 61 x 101 43 x 631 x 13567 43 x 271 x 5827 43 x 127 x 2731 43 x 127 x 1093 43 x 211 x 757 43 x 631 x 1597 43 x 127 x 211 43 x 211 x 337 43 x 433 x 643 43 x 547 x 673 43 x 3361 x 3907 47 x 3359 x 6073 47 x 1151 x 1933 47 x 3727 x 5153 53 x 157 x 2081 53 x 79 x 599 53 x 157 x 521 59 x 1451 x 2089 61 x 421 x 12841 61 x 181 x 5521 61 x 1301 x 19841 61 x 277 x 2113 61 x 181 x 1381 61 x 541 x 3001 61 x 661 x 2521 61 x 271 x 571 61 x 241 x 421 61 x 3361 x 4021
Haskell
<lang haskell>#!/usr/bin/runhaskell
import Data.Numbers.Primes import Data.Maybe
concatEach = flip concatMap eachMaybe = flip mapMaybe
carmichaels =
(takeWhile (<= 61) primes) `concatEach` (\p ->
[2..(p-1)] `concatEach` (\h3 ->
let g = h3 + p
in [1..(g-1)] `eachMaybe` (\d ->
if (g * (p - 1)) `mod` d == 0 && (-1 * p * p) `mod` h3 == d `mod` h3 then
let q = 1 + (((p - 1) * g) `div` d)
in if isPrime q then
let r = 1 + ((p * q) `div` h3)
in if isPrime r && (q * r) `mod` (p - 1) == 1 then
Just (p, q, r)
else
Nothing
else
Nothing
else
Nothing)))
main = putStr $ unlines $ map show carmichaels</lang>
- Output:
(3,11,17) (5,29,73) (5,17,29) (5,13,17) (7,19,67) (7,31,73) (7,13,31) (7,23,41) (7,73,103) (7,13,19) (13,61,397) (13,37,241) (13,97,421) (13,37,97) (13,37,61) (17,41,233) (17,353,1201) (19,43,409) (19,199,271) (23,199,353) (29,113,1093) (29,197,953) (31,991,15361) (31,61,631) (31,151,1171) (31,61,271) (31,61,211) (31,271,601) (31,181,331) (37,109,2017) (37,73,541) (37,613,1621) (37,73,181) (37,73,109) (41,1721,35281) (41,881,12041) (41,101,461) (41,241,761) (41,241,521) (41,73,137) (41,61,101) (43,631,13567) (43,271,5827) (43,127,2731) (43,127,1093) (43,211,757) (43,631,1597) (43,127,211) (43,211,337) (43,433,643) (43,547,673) (43,3361,3907) (47,3359,6073) (47,1151,1933) (47,3727,5153) (53,157,2081) (53,79,599) (53,157,521) (59,1451,2089) (61,421,12841) (61,181,5521) (61,1301,19841) (61,277,2113) (61,181,1381) (61,541,3001) (61,661,2521) (61,271,571) (61,241,421) (61,3361,4021)
Mathematica
<lang mathematica>Cases[Cases[
Cases[Table[{p1, h3, d}, {p1, Array[Prime, PrimePi@61]}, {h3, 2,
p1 - 1}, {d, 1, h3 + p1 - 1}], {p1_Integer, h3_, d_} /;
PrimeQ[1 + (p1 - 1) (h3 + p1)/d] &&
Divisible[p1^2 + d, h3] :> {p1, 1 + (p1 - 1) (h3 + p1)/d, h3},
Infinity], {p1_, p2_, h3_} /; PrimeQ[1 + Floor[p1 p2/h3]] :> {p1,
p2, 1 + Floor[p1 p2/h3]}], {p1_, p2_, p3_} /;
Mod[p2 p3, p1 - 1] == 1 :> Print[p1, "*", p2, "*", p3]]</lang>
Python
<lang python>class Isprime():
Extensible sieve of Erastosthenes
>>> isprime.check(11)
True
>>> isprime.multiples
{2, 4, 6, 8, 9, 10}
>>> isprime.primes
[2, 3, 5, 7, 11]
>>> isprime(13)
True
>>> isprime.multiples
{2, 4, 6, 7, 8, 9, 10, 11, 12, 14, 15, 16, 18, 20, 21, 22}
>>> isprime.primes
[2, 3, 5, 7, 11, 13, 17, 19]
>>> isprime.nmax
22
>>>
multiples = {2}
primes = [2]
nmax = 2
def __init__(self, nmax):
if nmax > self.nmax:
self.check(nmax)
def check(self, n):
if type(n) == float:
if not n.is_integer(): return False
n = int(n)
multiples = self.multiples
if n <= self.nmax:
return n not in multiples
else:
# Extend the sieve
primes, nmax = self.primes, self.nmax
newmax = max(nmax*2, n)
for p in primes:
multiples.update(range(p*((nmax + p + 1) // p), newmax+1, p))
for i in range(nmax+1, newmax+1):
if i not in multiples:
primes.append(i)
multiples.update(range(i*2, newmax+1, i))
self.nmax = newmax
return n not in multiples
__call__ = check
def carmichael(p1):
ans = []
if isprime(p1):
for h3 in range(2, p1):
g = h3 + p1
for d in range(1, g):
if (g * (p1 - 1)) % d == 0 and (-p1 * p1) % h3 == d % h3:
p2 = 1 + ((p1 - 1)* g // d)
if isprime(p2):
p3 = 1 + (p1 * p2 // h3)
if isprime(p3):
if (p2 * p3) % (p1 - 1) == 1:
#print('%i X %i X %i' % (p1, p2, p3))
ans += [tuple(sorted((p1, p2, p3)))]
return ans
isprime = Isprime(2)
ans = sorted(sum((carmichael(n) for n in range(62) if isprime(n)), [])) print(',\n'.join(repr(ans[i:i+5])[1:-1] for i in range(0, len(ans)+1, 5)))</lang>
- Output:
(3, 11, 17), (5, 13, 17), (5, 17, 29), (5, 29, 73), (7, 13, 19), (7, 13, 31), (7, 19, 67), (7, 23, 41), (7, 31, 73), (7, 73, 103), (13, 37, 61), (13, 37, 97), (13, 37, 241), (13, 61, 397), (13, 97, 421), (17, 41, 233), (17, 353, 1201), (19, 43, 409), (19, 199, 271), (23, 199, 353), (29, 113, 1093), (29, 197, 953), (31, 61, 211), (31, 61, 271), (31, 61, 631), (31, 151, 1171), (31, 181, 331), (31, 271, 601), (31, 991, 15361), (37, 73, 109), (37, 73, 181), (37, 73, 541), (37, 109, 2017), (37, 613, 1621), (41, 61, 101), (41, 73, 137), (41, 101, 461), (41, 241, 521), (41, 241, 761), (41, 881, 12041), (41, 1721, 35281), (43, 127, 211), (43, 127, 1093), (43, 127, 2731), (43, 211, 337), (43, 211, 757), (43, 271, 5827), (43, 433, 643), (43, 547, 673), (43, 631, 1597), (43, 631, 13567), (43, 3361, 3907), (47, 1151, 1933), (47, 3359, 6073), (47, 3727, 5153), (53, 79, 599), (53, 157, 521), (53, 157, 2081), (59, 1451, 2089), (61, 181, 1381), (61, 181, 5521), (61, 241, 421), (61, 271, 571), (61, 277, 2113), (61, 421, 12841), (61, 541, 3001), (61, 661, 2521), (61, 1301, 19841), (61, 3361, 4021)
REXX
Note that REXX's version of modulus (//) is really a remainder function, so a version of
the modulus function was hard-coded below (when using a negative value).
numbers in order of calculation
<lang rexx>/*REXX program calculates Carmichael 3-strong pseudoprimes (up to N).*/ numeric digits 30 /*in case user wants bigger nums.*/ parse arg N .; if N== then N=61 /*allow user to specify the limit*/ if 1=='f1'x then times='af'x /*if EBCDIC machine, use a bullet*/
else times='f9'x /* " ASCII " " " " */
carms=0 /*number of Carmichael #s so far.*/ !.=0 /*a method of prime memoization. */
do p=3 to N by 2; if \isPrime(p) then iterate /*Not prime? Skip.*/
pm=p-1; nps=-p*p; @.=0; min=1e9; max=0 /*some handy-dandy variables.*/
do h3=2 to pm; g=h3+p /*find Carmichael #s for this P. */
do d=1 to g-1
if g*pm//d\==0 then iterate
if ((nps//h3)+h3)//h3\==d//h3 then iterate
q=1+pm*g%d; if \isPrime(q) then iterate
r=1+p*q%h3; if q*r//pm\==1 then iterate
if \isPrime(r) then iterate
carms=carms+1 /*bump the Carmichael # counter. */
min=min(min,q); max=max(max,q); @.q=r /*build a list.*/
end /*d*/
end /*h3*/
/*display a list of some Carm #s.*/
do j=min to max by 2; if @.j==0 then iterate /*one of the #s?*/
say '──────── a Carmichael number: ' p times j times @.j
end /*j*/
say /*show bueatification blank line.*/
end /*p*/
say; say carms ' Carmichael numbers found.' exit /*stick a fork in it, we're done.*/ /*──────────────────────────────────ISPRIME subroutine──────────────────*/ isPrime: procedure expose !.; parse arg x; if !.x then return 1 if wordpos(x,'2 3 5 7 11 13')\==0 then do; !.x=1; return 1; end if x<17 then return 0; if x//2==0 then return 0; if x//3==0 then return 0 if right(x,1)==5 then return 0; if x//7==0 then return 0
do i=11 by 6 until i*i>x; if x// i ==0 then return 0
if x//(i+2) ==0 then return 0
end /*i*/
!.x=1; return 1</lang> output when using the default input
──────── a Carmichael number: 3 ∙ 11 ∙ 17 ──────── a Carmichael number: 5 ∙ 29 ∙ 73 ──────── a Carmichael number: 5 ∙ 17 ∙ 29 ──────── a Carmichael number: 5 ∙ 13 ∙ 17 ──────── a Carmichael number: 7 ∙ 19 ∙ 67 ──────── a Carmichael number: 7 ∙ 31 ∙ 73 ──────── a Carmichael number: 7 ∙ 13 ∙ 31 ──────── a Carmichael number: 7 ∙ 23 ∙ 41 ──────── a Carmichael number: 7 ∙ 73 ∙ 103 ──────── a Carmichael number: 7 ∙ 13 ∙ 19 ──────── a Carmichael number: 13 ∙ 61 ∙ 397 ──────── a Carmichael number: 13 ∙ 37 ∙ 241 ──────── a Carmichael number: 13 ∙ 97 ∙ 421 ──────── a Carmichael number: 13 ∙ 37 ∙ 97 ──────── a Carmichael number: 13 ∙ 37 ∙ 61 ──────── a Carmichael number: 17 ∙ 41 ∙ 233 ──────── a Carmichael number: 17 ∙ 353 ∙ 1201 ──────── a Carmichael number: 19 ∙ 43 ∙ 409 ──────── a Carmichael number: 19 ∙ 199 ∙ 271 ──────── a Carmichael number: 23 ∙ 199 ∙ 353 ──────── a Carmichael number: 29 ∙ 113 ∙ 1093 ──────── a Carmichael number: 29 ∙ 197 ∙ 953 ──────── a Carmichael number: 31 ∙ 991 ∙ 15361 ──────── a Carmichael number: 31 ∙ 61 ∙ 631 ──────── a Carmichael number: 31 ∙ 151 ∙ 1171 ──────── a Carmichael number: 31 ∙ 61 ∙ 271 ──────── a Carmichael number: 31 ∙ 61 ∙ 211 ──────── a Carmichael number: 31 ∙ 271 ∙ 601 ──────── a Carmichael number: 31 ∙ 181 ∙ 331 ──────── a Carmichael number: 37 ∙ 109 ∙ 2017 ──────── a Carmichael number: 37 ∙ 73 ∙ 541 ──────── a Carmichael number: 37 ∙ 613 ∙ 1621 ──────── a Carmichael number: 37 ∙ 73 ∙ 181 ──────── a Carmichael number: 37 ∙ 73 ∙ 109 ──────── a Carmichael number: 41 ∙ 1721 ∙ 35281 ──────── a Carmichael number: 41 ∙ 881 ∙ 12041 ──────── a Carmichael number: 41 ∙ 101 ∙ 461 ──────── a Carmichael number: 41 ∙ 241 ∙ 761 ──────── a Carmichael number: 41 ∙ 241 ∙ 521 ──────── a Carmichael number: 41 ∙ 73 ∙ 137 ──────── a Carmichael number: 41 ∙ 61 ∙ 101 ──────── a Carmichael number: 43 ∙ 631 ∙ 13567 ──────── a Carmichael number: 43 ∙ 271 ∙ 5827 ──────── a Carmichael number: 43 ∙ 127 ∙ 2731 ──────── a Carmichael number: 43 ∙ 127 ∙ 1093 ──────── a Carmichael number: 43 ∙ 211 ∙ 757 ──────── a Carmichael number: 43 ∙ 631 ∙ 1597 ──────── a Carmichael number: 43 ∙ 127 ∙ 211 ──────── a Carmichael number: 43 ∙ 211 ∙ 337 ──────── a Carmichael number: 43 ∙ 433 ∙ 643 ──────── a Carmichael number: 43 ∙ 547 ∙ 673 ──────── a Carmichael number: 43 ∙ 3361 ∙ 3907 ──────── a Carmichael number: 47 ∙ 3359 ∙ 6073 ──────── a Carmichael number: 47 ∙ 1151 ∙ 1933 ──────── a Carmichael number: 47 ∙ 3727 ∙ 5153 ──────── a Carmichael number: 53 ∙ 157 ∙ 2081 ──────── a Carmichael number: 53 ∙ 79 ∙ 599 ──────── a Carmichael number: 53 ∙ 157 ∙ 521 ──────── a Carmichael number: 59 ∙ 1451 ∙ 2089 ──────── a Carmichael number: 61 ∙ 421 ∙ 12841 ──────── a Carmichael number: 61 ∙ 181 ∙ 5521 ──────── a Carmichael number: 61 ∙ 1301 ∙ 19841 ──────── a Carmichael number: 61 ∙ 277 ∙ 2113 ──────── a Carmichael number: 61 ∙ 181 ∙ 1381 ──────── a Carmichael number: 61 ∙ 541 ∙ 3001 ──────── a Carmichael number: 61 ∙ 661 ∙ 2521 ──────── a Carmichael number: 61 ∙ 271 ∙ 571 ──────── a Carmichael number: 61 ∙ 241 ∙ 421 ──────── a Carmichael number: 61 ∙ 3361 ∙ 4021 69 Carmichael numbers found.
numbers in sorted order
With a few lines of code (and using sparse arrays), the Carmichael numbers can be shown in order. <lang rexx>/*REXX program calculates Carmichael 3-strong pseudoprimes (up to N).*/ numeric digits 30 /*in case user wants bigger nums.*/ parse arg N .; if N== then N=61 /*allow user to specify the limit*/ if 1=='f1'x then times='af'x /*if EBCDIC machine, use a bullet*/
else times='f9'x /* " ASCII " " " " */
carms=0 /*number of Carmichael #s so far.*/ !.=0 /*a method of prime memoization. */
/*Carmichael numbers aren't even.*/ do p=3 to N by 2; if \isPrime(p) then iterate /*Not prime? Skip.*/ pm=p-1; nps=-p*p; @.=0; min=1e9; max=0 /*some handy-dandy variables.*/
do h3=2 to pm; g=h3+p /*find Carmichael #s for this P. */
do d=1 to g-1
if g*pm//d\==0 then iterate
if ((nps//h3)+h3)//h3\==d//h3 then iterate
q=1+pm*g%d; if \isPrime(q) then iterate
r=1+p*q%h3; if q*r//pm\==1 then iterate
if \isPrime(r) then iterate
carms=carms+1 /*bump the Carmichael # counter. */
min=min(min,q); max=max(max,q); @.q=r /*build a list.*/
end /*d*/
end /*h3*/
/*display a list of some Carm #s.*/
do j=min to max by 2; if @.j==0 then iterate /*one of the #s?*/
say '──────── a Carmichael number: ' p times j times @.j
end /*j*/
say /*show bueatification blank line.*/
end /*p*/
say; say carms ' Carmichael numbers found.' exit /*stick a fork in it, we're done.*/ /*──────────────────────────────────ISPRIME subroutine──────────────────*/ isPrime: procedure expose !.; parse arg x; if !.x then return 1 if wordpos(x,'2 3 5 7 11 13')\==0 then do; !.x=1; return 1; end if x<17 then return 0; if x//2==0 then return 0; if x//3==0 then return 0 if right(x,1)==5 then return 0; if x//7==0 then return 0
do i=11 by 6 until i*i>x; if x// i ==0 then return 0
if x//(i+2) ==0 then return 0
end /*i*/
!.x=1; return 1</lang> output when using the default input
──────── a Carmichael number: 3 ∙ 11 ∙ 17 ──────── a Carmichael number: 5 ∙ 13 ∙ 17 ──────── a Carmichael number: 5 ∙ 17 ∙ 29 ──────── a Carmichael number: 5 ∙ 29 ∙ 73 ──────── a Carmichael number: 7 ∙ 13 ∙ 19 ──────── a Carmichael number: 7 ∙ 19 ∙ 67 ──────── a Carmichael number: 7 ∙ 23 ∙ 41 ──────── a Carmichael number: 7 ∙ 31 ∙ 73 ──────── a Carmichael number: 7 ∙ 73 ∙ 103 ──────── a Carmichael number: 13 ∙ 37 ∙ 61 ──────── a Carmichael number: 13 ∙ 61 ∙ 397 ──────── a Carmichael number: 13 ∙ 97 ∙ 421 ──────── a Carmichael number: 17 ∙ 41 ∙ 233 ──────── a Carmichael number: 17 ∙ 353 ∙ 1201 ──────── a Carmichael number: 19 ∙ 43 ∙ 409 ──────── a Carmichael number: 19 ∙ 199 ∙ 271 ──────── a Carmichael number: 23 ∙ 199 ∙ 353 ──────── a Carmichael number: 29 ∙ 113 ∙ 1093 ──────── a Carmichael number: 29 ∙ 197 ∙ 953 ──────── a Carmichael number: 31 ∙ 61 ∙ 211 ──────── a Carmichael number: 31 ∙ 151 ∙ 1171 ──────── a Carmichael number: 31 ∙ 181 ∙ 331 ──────── a Carmichael number: 31 ∙ 271 ∙ 601 ──────── a Carmichael number: 31 ∙ 991 ∙ 15361 ──────── a Carmichael number: 37 ∙ 73 ∙ 109 ──────── a Carmichael number: 37 ∙ 109 ∙ 2017 ──────── a Carmichael number: 37 ∙ 613 ∙ 1621 ──────── a Carmichael number: 41 ∙ 61 ∙ 101 ──────── a Carmichael number: 41 ∙ 73 ∙ 137 ──────── a Carmichael number: 41 ∙ 101 ∙ 461 ──────── a Carmichael number: 41 ∙ 241 ∙ 521 ──────── a Carmichael number: 41 ∙ 881 ∙ 12041 ──────── a Carmichael number: 41 ∙ 1721 ∙ 35281 ──────── a Carmichael number: 43 ∙ 127 ∙ 211 ──────── a Carmichael number: 43 ∙ 211 ∙ 337 ──────── a Carmichael number: 43 ∙ 271 ∙ 5827 ──────── a Carmichael number: 43 ∙ 433 ∙ 643 ──────── a Carmichael number: 43 ∙ 547 ∙ 673 ──────── a Carmichael number: 43 ∙ 631 ∙ 1597 ──────── a Carmichael number: 43 ∙ 3361 ∙ 3907 ──────── a Carmichael number: 47 ∙ 1151 ∙ 1933 ──────── a Carmichael number: 47 ∙ 3359 ∙ 6073 ──────── a Carmichael number: 47 ∙ 3727 ∙ 5153 ──────── a Carmichael number: 53 ∙ 79 ∙ 599 ──────── a Carmichael number: 53 ∙ 157 ∙ 521 ──────── a Carmichael number: 59 ∙ 1451 ∙ 2089 ──────── a Carmichael number: 61 ∙ 181 ∙ 1381 ──────── a Carmichael number: 61 ∙ 241 ∙ 421 ──────── a Carmichael number: 61 ∙ 271 ∙ 571 ──────── a Carmichael number: 61 ∙ 277 ∙ 2113 ──────── a Carmichael number: 61 ∙ 421 ∙ 12841 ──────── a Carmichael number: 61 ∙ 541 ∙ 3001 ──────── a Carmichael number: 61 ∙ 661 ∙ 2521 ──────── a Carmichael number: 61 ∙ 1301 ∙ 19841 ──────── a Carmichael number: 61 ∙ 3361 ∙ 4021 69 Carmichael numbers found.
Ruby
<lang ruby>
- Generate Charmichael Numbers
- Nigel_Galloway
- November 30th., 2012.
require 'prime' Integer.each_prime(61) {|p|
(2...p).each {|h3|
g = h3+p
(1...g).each {|d|
if (g*(p-1)) % d == 0 and (-1*p*p) % h3 == d % h3 then
q = 1+((p-1)*g/d)
next if not q.prime?
r = 1+(p*q/h3)
next if not r.prime? or not (q*r) % (p-1) == 1
puts "#{p} X #{q} X #{r}"
end
}
}
puts ""
} </lang>
- Output:
3 X 11 X 17 5 X 29 X 73 5 X 17 X 29 5 X 13 X 17 7 X 19 X 67 7 X 31 X 73 7 X 13 X 31 7 X 23 X 41 7 X 73 X 103 7 X 13 X 19 13 X 61 X 397 13 X 37 X 241 13 X 97 X 421 13 X 37 X 97 13 X 37 X 61 17 X 41 X 233 17 X 353 X 1201 19 X 43 X 409 19 X 199 X 271 23 X 199 X 353 29 X 113 X 1093 29 X 197 X 953 31 X 991 X 15361 31 X 61 X 631 31 X 151 X 1171 31 X 61 X 271 31 X 61 X 211 31 X 271 X 601 31 X 181 X 331 37 X 109 X 2017 37 X 73 X 541 37 X 613 X 1621 37 X 73 X 181 37 X 73 X 109 41 X 1721 X 35281 41 X 881 X 12041 41 X 101 X 461 41 X 241 X 761 41 X 241 X 521 41 X 73 X 137 41 X 61 X 101 43 X 631 X 13567 43 X 271 X 5827 43 X 127 X 2731 43 X 127 X 1093 43 X 211 X 757 43 X 631 X 1597 43 X 127 X 211 43 X 211 X 337 43 X 433 X 643 43 X 547 X 673 43 X 3361 X 3907 47 X 3359 X 6073 47 X 1151 X 1933 47 X 3727 X 5153 53 X 157 X 2081 53 X 79 X 599 53 X 157 X 521 59 X 1451 X 2089 61 X 421 X 12841 61 X 181 X 5521 61 X 1301 X 19841 61 X 277 X 2113 61 X 181 X 1381 61 X 541 X 3001 61 X 661 X 2521 61 X 271 X 571 61 X 241 X 421 61 X 3361 X 4021