CloudFlare suffered a massive security issue affecting all of its customers, including Rosetta Code. All passwords not changed since February 19th 2017 have been expired, and session cookie longevity will be reduced until late March.--Michael Mol (talk) 05:15, 25 February 2017 (UTC)

# Carmichael 3 strong pseudoprimes

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 separated from primes by Fermat's Little Theorem, but there are some that completely confound it.

The   Miller Rabin Test   uses a combination of Fermat's Little Theorem and Chinese Division Theorem to overcome this.

The purpose of this task is to investigate such numbers using a method based on   Carmichael numbers,   as suggested in   Notes by G.J.O Jameson March 2010.

Find Carmichael numbers of the form:

Prime1 × Prime2 × Prime3

where   (Prime1 < Prime2 < Prime3)   for all   Prime1   up to   61.
(See page 7 of   Notes by G.J.O Jameson March 2010   for solutions.)

Pseudocode

For a given   ${\displaystyle Prime_{1}}$

for 1 < 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) * (h3+Prime1)/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


Uses the Miller_Rabin package from Miller-Rabin primality test#ordinary integers.

with Ada.Text_IO, Miller_Rabin; procedure Nemesis is    type Number is range 0 .. 2**40-1; -- sufficiently large for the task    function Is_Prime(N: Number) return Boolean is      package MR is new Miller_Rabin(Number); use MR;   begin      return MR.Is_Prime(N) = Probably_Prime;   end Is_Prime; begin   for P1 in Number(2) .. 61 loop      if Is_Prime(P1) then         for H3 in Number(1) .. P1 loop            declare               G: Number := H3 + P1;               P2, P3: Number;            begin               Inner:               for D in 1 .. G-1 loop                  if ((H3+P1) * (P1-1)) mod D = 0 and then                    (-(P1 * P1)) mod H3 = D mod H3                  then                     P2 := 1 + ((P1-1) * G / D);                     P3 := 1 +(P1*P2/H3);                     if Is_Prime(P2) and then Is_Prime(P3)                       and then (P2*P3) mod (P1-1) = 1                     then                       Ada.Text_IO.Put_Line                        ( Number'Image(P1) & " *"   & Number'Image(P2) & " *" &                          Number'Image(P3) & "  = " & Number'Image(P1*P2*P3) );                     end if;                  end if;               end loop Inner;            end;         end loop;      end if;   end loop;end Nemesis;
Output:
 3 * 11 * 17  =  561
5 * 29 * 73  =  10585
5 * 17 * 29  =  2465
5 * 13 * 17  =  1105
7 * 19 * 67  =  8911

... (the full output is 69 lines long) ...

61 * 271 * 571  =  9439201
61 * 241 * 421  =  6189121
61 * 3361 * 4021  =  824389441

## ALGOL 68

Uses the Sieve of Eratosthenes code from the Smith Numbers task with an increased upper-bound (included here for convenience).

# sieve of Eratosthene: sets s[i] to TRUE if i is prime, FALSE otherwise #PROC sieve = ( REF[]BOOL s )VOID:     BEGIN        # start with everything flagged as prime                             #         FOR i TO UPB s DO s[ i ] := TRUE OD;        # sieve out the non-primes                                           #        s[ 1 ] := FALSE;        FOR i FROM 2 TO ENTIER sqrt( UPB s ) DO            IF s[ i ] THEN FOR p FROM i * i BY i TO UPB s DO s[ p ] := FALSE OD FI        OD     END # sieve # ; # construct a sieve of primes up to the maximum number required for the task ## For Prime1, we need to check numbers up to around 120 000                  #INT max number = 200 000;[ 1 : max number ]BOOL is prime;sieve( is prime ); # Find the Carmichael 3 Stromg Pseudoprimes for Prime1 up to 61              # FOR prime1 FROM 2 TO 61 DO    IF is prime[ prime 1 ] THEN        FOR h3 TO prime1 - 1 DO            FOR d TO ( h3 + prime1 ) - 1 DO                IF   ( h3 + prime1 ) * ( prime1 - 1 ) MOD d = 0                AND ( - ( prime1 * prime1 ) ) MOD h3 = d MOD h3                THEN                    INT prime2 = 1 + ( ( prime1 - 1 ) * ( h3 + prime1 ) OVER d );                    IF is prime[ prime2 ] THEN                        INT prime3 = 1 + ( prime1 * prime2 OVER h3 );                        IF is prime[ prime3 ] THEN                             IF ( prime2 * prime3 ) MOD ( prime1 - 1 ) = 1 THEN                                print( ( whole( prime1, 0 ), " ", whole( prime2, 0 ), " ", whole( prime3, 0 ), newline ) )                            FI                        FI                    FI                FI            OD        OD    FIOD
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
...
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


## C

 #include <stdio.h> /* C's % operator actually calculates the remainder of a / b so we need a * small adjustment so it works as expected for negative values */#define mod(n,m) ((((n) % (m)) + (m)) % (m)) int is_prime(unsigned int n){    if (n <= 3) {        return n > 1;    }    else if (!(n % 2) || !(n % 3)) {        return 0;    }    else {        unsigned int i;        for (i = 5; i*i <= n; i += 6)            if (!(n % i) || !(n % (i + 2)))                return 0;        return 1;    }} void carmichael3(int p1){    if (!is_prime(p1)) return;     int h3, d, p2, p3;    for (h3 = 1; h3 < p1; ++h3) {        for (d = 1; d < h3 + p1; ++d) {            if ((h3 + p1)*(p1 - 1) % d == 0 && mod(-p1 * p1, h3) == d % h3) {                p2 = 1 + ((p1 - 1) * (h3 + p1)/d);                if (!is_prime(p2)) continue;                p3 = 1 + (p1 * p2 / h3);                if (!is_prime(p3) || (p2 * p3) % (p1 - 1) != 1) continue;                printf("%d %d %d\n", p1, p2, p3);            }        }    }} int main(void){    int p1;    for (p1 = 2; p1 < 62; ++p1)        carmichael3(p1);    return 0;}
Output:
3 11 17
5 29 73
5 17 29
5 13 17
7 19 67
7 31 73
.
.
.
61 181 1381
61 541 3001
61 661 2521
61 271 571
61 241 421
61 3361 4021


## Clojure

 (ns example  (:gen-class)) (defn prime? [n]  " Prime number test (using Java) "  (.isProbablePrime (biginteger n) 16)) (defn carmichael [p1]  " Triplets of Carmichael primes, with first element prime p1 "  (if (prime? p1)    (into [] (for [h3 (range 2 p1)          :let [g (+ h3 p1)]          d (range 1 g)          :when (and (= (mod (* g (dec p1)) d) 0)                     (= (mod (- (* p1 p1)) h3) (mod d h3)))          :let [p2 (inc (quot (* (dec p1) g) d))]          :when (prime? p2)          :let [p3 (inc (quot (* p1 p2) h3))]          :when (prime? p3)          :when (= (mod (* p2 p3) (dec p1)) 1)]         [p1 p2 p3])))) ; Generate Result(def numbers (mapcat carmichael (range 2 62)))(println (count numbers) "Carmichael numbers found:")(doseq [t numbers]  (println (format "%5d x %5d x %5d = %10d" (first t) (second t) (last t) (apply * t))))
Output:
69 Carmichael numbers found
3 x    11 x    17 =        561
5 x    29 x    73 =      10585
5 x    17 x    29 =       2465
5 x    13 x    17 =       1105
7 x    19 x    67 =       8911
7 x    31 x    73 =      15841
7 x    13 x    31 =       2821
7 x    23 x    41 =       6601
7 x    73 x   103 =      52633
7 x    13 x    19 =       1729
13 x    61 x   397 =     314821
13 x    37 x   241 =     115921
13 x    97 x   421 =     530881
13 x    37 x    97 =      46657
13 x    37 x    61 =      29341
17 x    41 x   233 =     162401
17 x   353 x  1201 =    7207201
19 x    43 x   409 =     334153
19 x   199 x   271 =    1024651
23 x   199 x   353 =    1615681
29 x   113 x  1093 =    3581761
29 x   197 x   953 =    5444489
31 x   991 x 15361 =  471905281
31 x    61 x   631 =    1193221
31 x   151 x  1171 =    5481451
31 x    61 x   271 =     512461
31 x    61 x   211 =     399001
31 x   271 x   601 =    5049001
31 x   181 x   331 =    1857241
37 x   109 x  2017 =    8134561
37 x    73 x   541 =    1461241
37 x   613 x  1621 =   36765901
37 x    73 x   181 =     488881
37 x    73 x   109 =     294409
41 x  1721 x 35281 = 2489462641
41 x   881 x 12041 =  434932961
41 x   101 x   461 =    1909001
41 x   241 x   761 =    7519441
41 x   241 x   521 =    5148001
41 x    73 x   137 =     410041
41 x    61 x   101 =     252601
43 x   631 x 13567 =  368113411
43 x   271 x  5827 =   67902031
43 x   127 x  2731 =   14913991
43 x   127 x  1093 =    5968873
43 x   211 x   757 =    6868261
43 x   631 x  1597 =   43331401
43 x   127 x   211 =    1152271
43 x   211 x   337 =    3057601
43 x   433 x   643 =   11972017
43 x   547 x   673 =   15829633
43 x  3361 x  3907 =  564651361
47 x  3359 x  6073 =  958762729
47 x  1151 x  1933 =  104569501
47 x  3727 x  5153 =  902645857
53 x   157 x  2081 =   17316001
53 x    79 x   599 =    2508013
53 x   157 x   521 =    4335241
59 x  1451 x  2089 =  178837201
61 x   421 x 12841 =  329769721
61 x   181 x  5521 =   60957361
61 x  1301 x 19841 = 1574601601
61 x   277 x  2113 =   35703361
61 x   181 x  1381 =   15247621
61 x   541 x  3001 =   99036001
61 x   661 x  2521 =  101649241
61 x   271 x   571 =    9439201
61 x   241 x   421 =    6189121
61 x  3361 x  4021 =  824389441



## D

enum mod = (in int n, in int m) pure nothrow @nogc=> ((n % m) + m) % m; bool isPrime(in uint n) pure nothrow @nogc {  if (n == 2 || n == 3)    return true;  else if (n < 2 || n % 2 == 0 || n % 3 == 0)    return false;  for (uint div = 5, inc = 2; div ^^ 2 <= n;     div += inc, inc = 6 - inc)    if (n % div == 0)      return false;  return true;} void main() {  import std.stdio;   foreach (immutable p; 2 .. 62) {    if (!p.isPrime) continue;    foreach (immutable h3; 2 .. p) {      immutable g = h3 + p;      foreach (immutable d; 1 .. g) {        if ((g * (p - 1)) % d != 0 || mod(-p * p, h3) != d % h3)          continue;        immutable q = 1 + (p - 1) * g / d;        if (!q.isPrime) continue;        immutable r = 1 + (p * q / h3);        if (!r.isPrime || (q * r) % (p - 1) != 1) continue;        writeln(p, " x ", q, " x ", r);      }    }  }}
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

## EchoLisp

 ;; charmichaΓ«l numbers up to N-th prime ; 61 is 18-th prime(define (charms (N 18) local: (h31 0) (Prime2 0) (Prime3 0))(for* ((Prime1 (primes N))       (h3 (in-range 1 Prime1))       (d  (+ h3 Prime1)))      (set! h31 (+ h3 Prime1))      #:continue (!zero? (modulo (* h31 (1- Prime1)) d))      #:continue (!= (modulo d h3) (modulo (- (* Prime1 Prime1)) h3))      (set! Prime2 (1+ ( * (1- Prime1) (quotient h31 d))))      #:when (prime? Prime2)      (set! Prime3 (1+ (quotient (*  Prime1  Prime2)  h3)))      #:when (prime? Prime3)      #:when (= 1 (modulo (* Prime2 Prime3) (1- Prime1)))      (printf " π₯ %12d = %d x %d x %d"  (* Prime1 Prime2 Prime3) Prime1 Prime2 Prime3)))
Output:
 (charms 3)π₯          561 = 3 x 11 x 17π₯        10585 = 5 x 29 x 73π₯         2465 = 5 x 17 x 29π₯         1105 = 5 x 13 x 17 (charms 18);; skipped ....π₯    902645857 = 47 x 3727 x 5153π₯      2632033 = 53 x 53 x 937π₯     17316001 = 53 x 157 x 2081π₯      4335241 = 53 x 157 x 521π₯    178837201 = 59 x 1451 x 2089π₯    329769721 = 61 x 421 x 12841π₯     60957361 = 61 x 181 x 5521π₯      6924781 = 61 x 61 x 1861π₯      6924781 = 61 x 61 x 1861π₯     15247621 = 61 x 181 x 1381π₯     99036001 = 61 x 541 x 3001π₯    101649241 = 61 x 661 x 2521π₯      6189121 = 61 x 241 x 421π₯    824389441 = 61 x 3361 x 4021

## Fortran

### Plan

This is F77 style, and directly translates the given calculation as per formula translation. It turns out that the normal integers suffice for the demonstration, except for just one of the products of the three primes: 41x1721x35281 = 2489462641, which is bigger than 2147483647, the 32-bit limit. Fortunately, INTEGER*8 variables are also available, so the extension is easy. Otherwise, one would have to mess about with using two integers in a bignum style, one holding say the millions, and the second the number up to a million.

### Integer arithmetic

Integer arithmetic can be a delicate matter, not just because of overflows. With division, sometimes the formula will never produce a remainder, as in n(n + 1)/2, and sometimes the divisor is always a factor for more complex reasons - as with (P1 - 1)*(H3 + P1)/D. But if a term produces a fractional part, is it to be discarded or not? Some languages (such as Pascal) insist on "div" instead of "/" so that one thinks about truncating division, others (such as pl/i) just generate the fractional part and proceed: perhaps the result will round? truncate? to the desired value, or perhaps other terms will cancel out the fractional parts and all will be well, or, perhaps they will combine and the result will be out by one, etc. And some formulae require the fractional parts, as in the expression for the n'th Fibonacci number where each term generates an infinite non-repeating fractional part, and, they cancel exactly, even for large values of n... In source code,

(((1 + SQRT(5))/2)**N - ((1 - SQRT(5))/2)**N)/SQRT(5)

or, in mathematical notation (possibly not rendered),

${\displaystyle {\frac {1}{\sqrt {5}}}\left[\left({\frac {1+{\sqrt {5}}}{2}}\right)^{n}-\left({\frac {1-{\sqrt {5}}}{2}}\right)^{n}\right]}$

### What model MOD is your MOD?

Directly translating the given formulae involves a trap. There is an expression -Prime1 squared mod h3, and translating this to MOD(-P1**2,H3) may not work. The behaviour of the MOD function for negative numbers varies between language, compiler, and cpu. For positive numbers there is no confusion. However, in some cases one wants an affine shift as in day-of-the-week calculations from a daynumber, D: given MOD(D,7), one wants the same result with MOD(D - 7,7) or any other such shift and indeed, this happens with some versions of MOD. But not all, as Prof. Knuth remarks in his description of the calculation of the date for Easter. The MOD function can also work as follows:

             D        -7 -6 -5 -4 -3 -2 -1  0  1  2  3  4  5  6  7
MOD(D,3)     -1  0 -2 -1  0 -2 -1  0  1  2  0  1  2  0  1
MOD(3 + MOD(D,3),3)   2  0  1  2  0  1  2  0  1  2  0  1  2  0  1


Thus, MOD(-1,3) + MOD(4,3) = -1 + 1 and MOD(3,3) = 0 so that (a + b) mod n does come out as [(a mod n) + (b mod n)] mod n as is expected, but in a different way: (-1 mod 3) + (4 mod 3) = 2 + 1, thus the need for the second mod. So, MOD and mod are different. However, this MOD is consistent with integer division that truncates. If qΒ·f = n/d where q is the integer part and f the fractional part of the division, then the remainder is fd and has a sign too. Another computer/compiler/language version of MOD may not differ from mod.

Here I am supposing that (h3+Prime1)*(Prime1-1) mod d == 0 and -Prime1 squared mod h3 == d mod h3 is evaluated as {[(h3 + Prime1)*(Prime1 - 1) mod d] == 0} and {[-(Prime1 squared) mod h3] == d mod h3} which is to say it is [-(Prime1 squared)] mod h3 not -[(Prime1 squared) mod h3] where the first case involves the MOD(negative number) while the second is -MOD(positive number) and the only way that result could match d mod h3 is when both are zero.

Testing with the peculiar sign ignored via .AND. (MOD(P1**2,H3) .EQ. MOD(D,H3))) THEN gives limited results:

Carmichael numbers that are the product of three primes:
P1  x P2  x P3 =         C
3    11    17         561
5    29    73       10585
7    19    67        8911
13    61   397      314821
13    37   241      115921
19    43   409      334153
31   991 15361   471905281
37   109  2017     8134561
41  1721 35281  2489462641
43   631 13567   368113411
43   271  5827    67902031
43   127  2731    14913991
61   421 12841   329769721
61   181  5521    60957361

Which turns out to be when the two terms have remainders of one.

### Source

So, using the double MOD approach - which gives the same result for either style of MOD...
      LOGICAL FUNCTION ISPRIME(N)	!Ad-hoc, since N is not going to be big...       INTEGER N			!Despite this intimidating allowance of 32 bits...       INTEGER F			!A possible factor.        ISPRIME = .FALSE.		!Most numbers aren't prime.        DO F = 2,SQRT(DFLOAT(N))	!Wince...          IF (MOD(N,F).EQ.0) RETURN	!Not even avoiding even numbers beyond two.        END DO				!Nice and brief, though.        ISPRIME = .TRUE.		!No factor found.      END FUNCTION ISPRIME		!So, done. Hopefully, not often.       PROGRAM CHASE      INTEGER P1,P2,P3	!The three primes to be tested.      INTEGER H3,D	!Assistants.      INTEGER MSG	!File unit number.      MSG = 6		!Standard output.      WRITE (MSG,1)	!A heading would be good.    1 FORMAT ("Carmichael numbers that are the product of three primes:"     & /"    P1  x P2  x P3 =",9X,"C")      DO P1 = 2,61	!Step through the specified range.        IF (ISPRIME(P1)) THEN	!Selecting only the primes.          DO H3 = 2,P1 - 1		!For 1 < H3 < P1.            DO D = 1,H3 + P1 - 1		!For 0 < D < H3 + P1.              IF (MOD((H3 + P1)*(P1 - 1),D).EQ.0	!Filter.     &        .AND. (MOD(H3 + MOD(-P1**2,H3),H3) .EQ. MOD(D,H3))) THEN	!Beware MOD for negative numbers! MOD(-P1**2, may surprise...                P2 = 1 + (P1 - 1)*(H3 + P1)/D	!Candidate for the second prime.                IF (ISPRIME(P2)) THEN		!Is it prime?                  P3 = 1 + P1*P2/H3			!Yes. Candidate for the third prime.                  IF (ISPRIME(P3)) THEN			!Is it prime?                    IF (MOD(P2*P3,P1 - 1).EQ.1) THEN		!Yes! Final test.                      WRITE (MSG,2) P1,P2,P3, INT8(P1)*P2*P3		!Result!    2                 FORMAT (3I6,I12)                    END IF                  END IF                END IF              END IF            END DO          END DO        END IF      END DO      END

### Output

Carmichael numbers that are the product of three primes:
P1  x P2  x P3 =         C
3    11    17         561
5    29    73       10585
5    17    29        2465
5    13    17        1105
7    19    67        8911
7    31    73       15841
7    13    31        2821
7    23    41        6601
7    73   103       52633
7    13    19        1729
13    61   397      314821
13    37   241      115921
13    97   421      530881
13    37    97       46657
13    37    61       29341
17    41   233      162401
17   353  1201     7207201
19    43   409      334153
19   199   271     1024651
23   199   353     1615681
29   113  1093     3581761
29   197   953     5444489
31   991 15361   471905281
31    61   631     1193221
31   151  1171     5481451
31    61   271      512461
31    61   211      399001
31   271   601     5049001
31   181   331     1857241
37   109  2017     8134561
37    73   541     1461241
37   613  1621    36765901
37    73   181      488881
37    73   109      294409
41  1721 35281  2489462641
41   881 12041   434932961
41   101   461     1909001
41   241   761     7519441
41   241   521     5148001
41    73   137      410041
41    61   101      252601
43   631 13567   368113411
43   271  5827    67902031
43   127  2731    14913991
43   127  1093     5968873
43   211   757     6868261
43   631  1597    43331401
43   127   211     1152271
43   211   337     3057601
43   433   643    11972017
43   547   673    15829633
43  3361  3907   564651361
47  3359  6073   958762729
47  1151  1933   104569501
47  3727  5153   902645857
53   157  2081    17316001
53    79   599     2508013
53   157   521     4335241
59  1451  2089   178837201
61   421 12841   329769721
61   181  5521    60957361
61  1301 19841  1574601601
61   277  2113    35703361
61   181  1381    15247621
61   541  3001    99036001
61   661  2521   101649241
61   271   571     9439201
61   241   421     6189121
61  3361  4021   824389441


## FreeBASIC

' version 17-10-2016' compile with: fbc -s console ' using a sieve for finding primes #Define max_sieve 10000000 ' 10^7ReDim Shared As Byte isprime(max_sieve) ' translated the pseudo code to FreeBASIC Sub carmichael3(p1 As Integer)    If isprime(p1) = 0 Then Exit Sub   Dim As Integer h3, d, p2, p3, t1, t2   For h3 = 1 To p1 -1    t1 = (h3 + p1) * (p1 -1)    t2 = (-p1 * p1) Mod h3    If t2 < 0 Then t2 = t2 + h3    For d = 1 To h3 + p1 -1      If t1 Mod d = 0 And t2 = (d Mod h3) Then        p2 = 1 + (t1 \ d)        If isprime(p2) = 0 Then Continue For        p3 = 1 + (p1 * p2 \ h3)        If isprime(p3) = 0 Or ((p2 * p3) Mod (p1 -1)) <> 1 Then Continue For        Print Using "### * #### * #####"; p1; p2; p3      End If    Next d  Next h3End Sub  ' ------=< MAIN >=------ Dim As UInteger i, j 'set up sieveFor i = 3 To max_sieve Step 2  isprime(i) = 1Next i isprime(2) = 1For i = 3 To Sqr(max_sieve) Step 2  If isprime(i) = 1 Then    For j = i * i To max_sieve Step i * 2      isprime(j) = 0    Next j  End IfNext i For i = 2 To 61  carmichael3(i)Next i ' empty keyboard bufferWhile InKey <> "" : WendPrint : Print "hit any key to end program"SleepEnd
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

Translation of: Ruby
Library: primes
Works with: GHC version 7.4.1
Works with: primes version 0.2.1.0
Output:

See D output.

## Mathematica / Wolfram Language

Cases[Cases[  Cases[Table[{p1, h3, d}, {p1, Array[Prime, [email protected]/* <![CDATA[ */!function(t,e,r,n,c,a,p){try{t=document.currentScript||function(){for(t=document.getElementsByTagName('script'),e=t.length;e--;)if(t[e].getAttribute('data-cfhash'))return t[e]}();if(t&&(c=t.previousSibling)){p=t.parentNode;if(a=c.getAttribute('data-cfemail')){for(e='',r='0x'+a.substr(0,2)|0,n=2;a.length-n;n+=2)e+='%'+('0'+('0x'+a.substr(n,2)^r).toString(16)).slice(-2);p.replaceChild(document.createTextNode(decodeURIComponent(e)),c)}p.removeChild(t)}}catch(u){}}()/* ]]> */]}, {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]]

## PARI/GP

f(p)={  my(v=List(),q,r);  for(h=2,p-1,    for(d=1,h+p-1,      if((h+p)*(p-1)%d==0 && Mod(p,h)^2==-d && isprime(q=(p-1)*(h+p)/d+1) && isprime(r=p*q\h+1)&&q*r%(p-1)==1,        listput(v,p*q*r)      )    )  );  Set(v)};forprime(p=3,67,v=f(p); for(i=1,#v,print1(v[i]", ")))
Output:
561, 1105, 2465, 10585, 1729, 2821, 6601, 8911, 15841, 52633, 29341, 46657, 115921, 314821, 530881, 162401, 7207201, 334153, 1024651, 1615681, 3581761, 5444489, 399001, 512461, 1193221, 1857241, 5049001, 5481451, 471905281, 294409, 488881, 1461241, 8134561, 36765901, 252601, 410041, 1909001, 5148001, 7519441, 434932961, 2489462641, 1152271, 3057601, 5968873, 6868261, 11972017, 14913991, 15829633, 43331401, 67902031, 368113411, 564651361, 104569501, 902645857, 958762729, 2508013, 4335241, 17316001, 178837201, 6189121, 9439201, 15247621, 35703361, 60957361, 99036001, 101649241, 329769721, 824389441, 1574601601, 10267951, 163954561, 7991602081,

## Perl

Library: ntheory
use ntheory qw/forprimes is_prime vecprod/; forprimes { my $p =$_;   for my $h3 (2 ..$p-1) {      my $ph3 =$p + $h3; for my$d (1 .. $ph3-1) { # Jameseon procedure page 6 next if ((-$p*$p) %$h3) != ($d %$h3);         next if (($p-1)*$ph3) % $d; my$q = 1 + ($p-1)*$ph3 / $d; # Jameson eq 7 next unless is_prime($q);         my $r = 1 + ($p*$q-1) /$h3;         # Jameson eq 6         next unless is_prime($r); next unless ($q*$r) % ($p-1) == 1;         printf "%2d x %5d x %5d = %s\n",$p,$q,$r,vecprod($p,$q,$r);      }   }} 3,61;
Output:
 3 x    11 x    17 = 561
5 x    29 x    73 = 10585
5 x    17 x    29 = 2465
5 x    13 x    17 = 1105
... full output is 69 lines ...
61 x   661 x  2521 = 101649241
61 x   271 x   571 = 9439201
61 x   241 x   421 = 6189121
61 x  3361 x  4021 = 824389441


## Perl 6

Works with: Rakudo version 2015.12

An almost direct translation of the pseudocode. We take the liberty of going up to 67 to show we aren't limited to 32-bit integers. (Perl 6 uses arbitrary precision in any case.)

for (2..67).grep: *.is-prime -> \Prime1 {    for 1 ^..^ Prime1 -> \h3 {        my \g = h3 + Prime1;        for 0 ^..^ h3 + Prime1 -> \d {            if (h3 + Prime1) * (Prime1 - 1) %% d and -Prime1**2 % h3 == d % h3  {                my \Prime2 = floor 1 + (Prime1 - 1) * g / d;                next unless Prime2.is-prime;                my \Prime3 = floor 1 + Prime1 * Prime2 / h3;                next unless Prime3.is-prime;                next unless (Prime2 * Prime3) % (Prime1 - 1) == 1;                say "{Prime1} Γ {Prime2} Γ {Prime3} == {Prime1 * Prime2 * Prime3}";            }        }    }}
Output:
3 Γ 11 Γ 17 == 561
5 Γ 29 Γ 73 == 10585
5 Γ 17 Γ 29 == 2465
5 Γ 13 Γ 17 == 1105
7 Γ 19 Γ 67 == 8911
7 Γ 31 Γ 73 == 15841
7 Γ 13 Γ 31 == 2821
7 Γ 23 Γ 41 == 6601
7 Γ 73 Γ 103 == 52633
7 Γ 13 Γ 19 == 1729
13 Γ 61 Γ 397 == 314821
13 Γ 37 Γ 241 == 115921
13 Γ 97 Γ 421 == 530881
13 Γ 37 Γ 97 == 46657
13 Γ 37 Γ 61 == 29341
17 Γ 41 Γ 233 == 162401
17 Γ 353 Γ 1201 == 7207201
19 Γ 43 Γ 409 == 334153
19 Γ 199 Γ 271 == 1024651
23 Γ 199 Γ 353 == 1615681
29 Γ 113 Γ 1093 == 3581761
29 Γ 197 Γ 953 == 5444489
31 Γ 991 Γ 15361 == 471905281
31 Γ 61 Γ 631 == 1193221
31 Γ 151 Γ 1171 == 5481451
31 Γ 61 Γ 271 == 512461
31 Γ 61 Γ 211 == 399001
31 Γ 271 Γ 601 == 5049001
31 Γ 181 Γ 331 == 1857241
37 Γ 109 Γ 2017 == 8134561
37 Γ 73 Γ 541 == 1461241
37 Γ 613 Γ 1621 == 36765901
37 Γ 73 Γ 181 == 488881
37 Γ 73 Γ 109 == 294409
41 Γ 1721 Γ 35281 == 2489462641
41 Γ 881 Γ 12041 == 434932961
41 Γ 101 Γ 461 == 1909001
41 Γ 241 Γ 761 == 7519441
41 Γ 241 Γ 521 == 5148001
41 Γ 73 Γ 137 == 410041
41 Γ 61 Γ 101 == 252601
43 Γ 631 Γ 13567 == 368113411
43 Γ 271 Γ 5827 == 67902031
43 Γ 127 Γ 2731 == 14913991
43 Γ 127 Γ 1093 == 5968873
43 Γ 211 Γ 757 == 6868261
43 Γ 631 Γ 1597 == 43331401
43 Γ 127 Γ 211 == 1152271
43 Γ 211 Γ 337 == 3057601
43 Γ 433 Γ 643 == 11972017
43 Γ 547 Γ 673 == 15829633
43 Γ 3361 Γ 3907 == 564651361
47 Γ 3359 Γ 6073 == 958762729
47 Γ 1151 Γ 1933 == 104569501
47 Γ 3727 Γ 5153 == 902645857
53 Γ 157 Γ 2081 == 17316001
53 Γ 79 Γ 599 == 2508013
53 Γ 157 Γ 521 == 4335241
59 Γ 1451 Γ 2089 == 178837201
61 Γ 421 Γ 12841 == 329769721
61 Γ 181 Γ 5521 == 60957361
61 Γ 1301 Γ 19841 == 1574601601
61 Γ 277 Γ 2113 == 35703361
61 Γ 181 Γ 1381 == 15247621
61 Γ 541 Γ 3001 == 99036001
61 Γ 661 Γ 2521 == 101649241
61 Γ 271 Γ 571 == 9439201
61 Γ 241 Γ 421 == 6189121
61 Γ 3361 Γ 4021 == 824389441
67 Γ 2311 Γ 51613 == 7991602081
67 Γ 331 Γ 7393 == 163954561
67 Γ 331 Γ 463 == 10267951

## PicoLisp

(de modulo (X Y)   (% (+ Y (% X Y)) Y) ) (de prime? (N)   (let D 0      (or         (= N 2)         (and            (> N 1)            (bit? 1 N)            (for (D 3  T  (+ D 2))               (T (> D (sqrt N)) T)               (T (=0 (% N D)) NIL) ) ) ) ) ) (for P1 61   (when (prime? P1)      (for (H3 2 (> P1 H3) (inc H3))         (let G (+ H3 P1)            (for (D 1 (> G D) (inc D))               (when                  (and                     (=0                        (% (* G (dec P1)) D) )                     (=                        (modulo (* (- P1) P1) H3)                        (% D H3)) )                  (let                     (P2                        (inc                           (/ (* (dec P1) G) D) )                        P3 (inc (/ (* P1 P2) H3)) )                     (when                        (and                           (prime? P2)                           (prime? P3)                           (= 1 (modulo (* P2 P3) (dec P1))) )                        (print (list P1 P2 P3)) ) ) ) ) ) ) ) )(prinl) (bye)

## PL/I

Carmichael: procedure options (main, reorder);  /* 24 January 2014 */   declare (Prime1, Prime2, Prime3, h3, d) fixed binary (31);    put ('Carmichael numbers are:');    do Prime1 = 1 to 61;       do h3 = 2 to Prime1; d_loop:  do d = 1 to h3+Prime1-1;            if (mod((h3+Prime1)*(Prime1-1), d) = 0) &               (mod(-Prime1*Prime1, h3) = mod(d, h3)) then               do;                  Prime2 = (Prime1-1) * (h3+Prime1)/d; Prime2 = Prime2 + 1;                  if ^is_prime(Prime2) then iterate d_loop;                  Prime3 = Prime1*Prime2/h3; Prime3 = Prime3 + 1;                  if ^is_prime(Prime3) then iterate d_loop;                  if mod(Prime2*Prime3, Prime1-1) ^= 1 then iterate d_loop;                  put skip edit (trim(Prime1), ' x ', trim(Prime2), ' x ', trim(Prime3)) (A);               end;         end;      end;   end;    /* Uses is_prime from Rosetta Code PL/I. */ end Carmichael;

Results:

Carmichael numbers are:
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
9 x 89 x 401
9 x 29 x 53
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
21 x 761 x 941
23 x 199 x 353
27 x 131 x 443
27 x 53 x 131
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
35 x 647 x 7549
35 x 443 x 3877
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
49 x 313 x 5113
49 x 97 x 433
51 x 701 x 7151
53 x 157 x 2081
53 x 79 x 599
53 x 157 x 521
55 x 3079 x 84673
55 x 163 x 4483
55 x 1567 x 28729
55 x 109 x 1999
55 x 433 x 2647
55 x 919 x 3889
55 x 139 x 547
55 x 3889 x 12583
55 x 109 x 163
55 x 433 x 487
57 x 113 x 1289
57 x 113 x 281
57 x 4649 x 10193
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


## Python

class Isprime():    '''    Extensible sieve of Eratosthenes     >>> 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)))
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)

## Racket

 #lang racket(require math) (for ([p1 (in-range 3 62)] #:when (prime? p1))  (for ([h3 (in-range 2 p1)])    (define g (+ p1 h3))    (let next ([d 1])      (when (< d g)        (when (and (zero? (modulo (* g (- p1 1)) d))                   (= (modulo (- (sqr p1)) h3) (modulo d h3)))          (define p2 (+ 1 (quotient (* g (- p1 1)) d)))          (when (prime? p2)            (define p3 (+ 1 (quotient (* p1 p2) h3)))            (when (and (prime? p3) (= 1 (modulo (* p2 p3) (- p1 1))))              (displayln (list p1 p2 p3 '=> (* p1 p2 p3))))))        (next (+ d 1))))))

Output:

 (3 11 17 => 561)(5 29 73 => 10585)(5 17 29 => 2465)(5 13 17 => 1105)(7 19 67 => 8911)(7 31 73 => 15841)(7 23 41 => 6601)(7 73 103 => 52633)(13 61 397 => 314821)(13 97 421 => 530881)(13 37 97 => 46657)(13 37 61 => 29341)(17 41 233 => 162401)(17 353 1201 => 7207201)(19 43 409 => 334153)(19 199 271 => 1024651)(23 199 353 => 1615681)(29 113 1093 => 3581761)(29 197 953 => 5444489)(31 991 15361 => 471905281)(31 61 631 => 1193221)(31 151 1171 => 5481451)(31 61 271 => 512461)(31 61 211 => 399001)(31 271 601 => 5049001)(31 181 331 => 1857241)(37 109 2017 => 8134561)(37 73 541 => 1461241)(37 613 1621 => 36765901)(37 73 181 => 488881)(37 73 109 => 294409)(41 1721 35281 => 2489462641)(41 881 12041 => 434932961)(41 101 461 => 1909001)(41 241 761 => 7519441)(41 241 521 => 5148001)(41 73 137 => 410041)(41 61 101 => 252601)(43 631 13567 => 368113411)(43 127 1093 => 5968873)(43 211 757 => 6868261)(43 631 1597 => 43331401)(43 127 211 => 1152271)(43 211 337 => 3057601)(43 433 643 => 11972017)(43 547 673 => 15829633)(43 3361 3907 => 564651361)(47 3359 6073 => 958762729)(47 1151 1933 => 104569501)(47 3727 5153 => 902645857)(53 157 2081 => 17316001)(53 79 599 => 2508013)(53 157 521 => 4335241)(59 1451 2089 => 178837201)(61 421 12841 => 329769721)(61 1301 19841 => 1574601601)(61 277 2113 => 35703361)(61 541 3001 => 99036001)(61 661 2521 => 101649241)(61 271 571 => 9439201)(61 241 421 => 6189121)(61 3361 4021 => 824389441)

## REXX

### slightly optimized

Note that REXX's version of   modulus   (//)   is really a   remainder   function.

The Carmichael numbers are shown in numerical order.

Some code optimization was done, while not necessary for the small default number (61),   it was significant for larger numbers.

/*REXX program calculates  Carmichael  3βstrong  pseudoprimes  (up to and including N). */numeric digits 18                                /*handle big dig #s (9 is the default).*/parse arg N .;    if N==''  then N=61            /*allow user to specify for the search.*/tell= N>0;           N=abs(N)                    /*N>0?  Then display Carmichael numbers*/#=0                                              /*number of Carmichael numbers so far. */@.=0;   @.2=1; @.3=1; @.5=1; @.7=1; @.11=1; @.13=1; @.17=1; @.19=1; @.23=1; @.29=1; @.31=1                                                 /*[β]  prime number memoization array. */    do p=3  to N  by 2;  pm=p-1;   bot=0;  top=0 /*step through some (odd) prime numbers*/    if \isPrime(p)  then iterate;  nps=-p*p      /*is   P   a prime?   No, then skip it.*/    !.=0                                         /*the list of Carmichael #'s  (so far).*/             do h3=2  to  pm;  g=h3+p            /*find Carmichael #s  for this prime.  */             gPM=g*pm;  npsH3=((nps//h3)+h3)//h3 /*define a couple of shortcuts for pgm.*/                                                 /* [β] perform some weeding of D values*/                 do d=1  for g-1;                   if gPM//d \== 0      then iterate                                                    if npsH3  \== d//h3  then iterate                                       q=1+gPM%d;   if \isPrime(q)       then iterate                                       r=1+p*q%h3;  if q*r//pm\==1       then iterate                                                    if \isPrime(r)       then iterate                 #=#+1;                !.q=r     /*bump Carmichael counter; add to array*/                 if bot==0  then bot=q;   bot=min(bot,q);    top=max(top,q)                 end   /*d*/             end       /*h3*/    $= /*display a list of some Carmichael #s.*/ do j=bot to top by 2 while tell; if !.j\==0 then$=$p"β"j'β'!.j end /*j*/ if$\==''  then say 'Carmichael number: '      strip($) end /*p*/saysay 'ββββββββ ' # " Carmichael numbers found."exit /*stick a fork in it, we're all done. *//*ββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββ*/isPrime: parse arg x; if @.x then return 1 /*X a known prime?*/ if x<37 then return 0; if x//2==0 then return 0; if x// 3==0 then return 0 parse var x '' -1 _; if _==5 then return 0; if x// 7==0 then return 0 do k=11 by 6 until k*k>x; if x// k ==0 then return 0 if x//(k+2)==0 then return 0 end /*i*/ @.x=1; return 1 output when using the default input: Carmichael number: 3β11β17 Carmichael number: 5β13β17 5β17β29 5β29β73 Carmichael number: 7β13β19 7β19β67 7β23β41 7β31β73 7β73β103 Carmichael number: 13β37β61 13β61β397 13β97β421 Carmichael number: 17β41β233 17β353β1201 Carmichael number: 19β43β409 19β199β271 Carmichael number: 23β199β353 Carmichael number: 29β113β1093 29β197β953 Carmichael number: 31β61β211 31β151β1171 31β181β331 31β271β601 31β991β15361 Carmichael number: 37β73β109 37β109β2017 37β613β1621 Carmichael number: 41β61β101 41β73β137 41β101β461 41β241β521 41β881β12041 41β1721β35281 Carmichael number: 43β127β211 43β211β337 43β271β5827 43β433β643 43β547β673 43β631β1597 43β3361β3907 Carmichael number: 47β1151β1933 47β3359β6073 47β3727β5153 Carmichael number: 53β79β599 53β157β521 Carmichael number: 59β1451β2089 Carmichael number: 61β181β1381 61β241β421 61β271β571 61β277β2113 61β421β12841 61β541β3001 61β661β2521 61β1301β19841 61β3361β4021 ββββββββ 69 Carmichael numbers found.  output when using the input of: -1000 ββββββββ 1038 Carmichael numbers found.  output when using the input of: -10000 ββββββββ 8716 Carmichael numbers found.  ### more optimized This REXX version (pre-)generates a number of primes to assist the isPrime function. /*REXX program calculates Carmichael 3βstrong pseudoprimes (up to and including N). */numeric digits 18 /*handle big dig #s (9 is the default).*/parse arg N .; if N=='' then N=61 /*allow user to specify for the search.*/tell= N>0; N=abs(N) /*N>0? Then display Carmichael numbers*/#=0; @.=0 /*number of Carmichael numbers so far. */@.2=1; @.3=1; @.5=1; @.7=1; @.11=1; @.13=1; @.17=1; @.19=1; @.23=1; @.29=1; @.31=1; @.37=1HP=37; do i=HP+2 by 2 for N*20; if isPrime(i) then do; @.i=1; HP=i; end; end /*i*/HP=HP+2 /*[β] prime number memoization array. */ do p=3 to N by 2; pm=p-1; bot=0; top=0 /*step through some (odd) prime numbers*/ if \isPrime(p) then iterate; nps=-p*p /*is P a prime? No, then skip it.*/ !.=0 /*the list of Carmichael #'s (so far).*/ do h3=2 to pm; g=h3+p /*find Carmichael #s for this prime. */ gPM=g*pm; npsH3=((nps//h3)+h3)//h3 /*define a couple of shortcuts for pgm.*/ /* [β] perform some weeding of D values*/ do d=1 for g-1; if gPM//d \== 0 then iterate if npsH3 \== d//h3 then iterate q=1+gPM%d; if \isPrime(q) then iterate r=1+p*q%h3; if q*r//pm\==1 then iterate if \isPrime(r) then iterate #=#+1; !.q=r /*bump Carmichael counter; add to array*/ if bot==0 then bot=q; bot=min(bot,q); top=max(top,q) end /*d*/ end /*h3*/$=                                           /*display a list of some Carmichael #s.*/             do j=bot  to top  by 2  while tell;   if !.j\==0  then $=$  p"β"j'β'!.j             end           /*j*/     if $\=='' then say 'Carmichael number: ' strip($)    end                /*p*/saysay 'ββββββββ '     #     " Carmichael numbers found."exit                                             /*stick a fork in it,  we're all done. *//*ββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββββ*/isPrime: parse arg x;             if @.x      then return 1           /*X a known prime?*/         if x<HP  then return 0;  if x//2==0  then return 0; if x// 3==0     then return 0         parse var x '' -1 _;     if _==5     then return 0; if x// 7==0     then return 0         if x//11==0      then return 0;                     if x//13==0     then return 0         if x//17==0      then return 0;                     if x//19==0     then return 0         if x//23==0      then return 0;                     if x//29==0     then return 0         if x//31==0      then return 0;                     if x//37==0     then return 0c=x                                              /*                                  ___*/b=1;   do  while b<=x;  b=b*4;  end              /*these two lines compute integer  β X */s=0;   do  while b>1;   b=b%4;  _=c-s-b;  s=s%2;   if _>=0  then do; c=_; s=s+b; end;  end                                do k=41  by 6  to s;               parse var  k  ''  -1  _                                if _\==5  then               if x// k   ==0  then return 0                                if _\==3  then               if x//(k+2)==0  then return 0                                end  /*k*/       /*K   will never be divisible by three.*/         @.x=1;   return 1                       /*Define a new prime (X).  Indicate so.*/

output   when the following us used for input:   1000

Carmichael number:  3β11β17
Carmichael number:  5β13β17 5β17β29 5β29β73
Carmichael number:  7β13β19 7β19β67 7β23β41 7β31β73 7β73β103
Carmichael number:  13β37β61 13β61β397 13β97β421
Carmichael number:  17β41β233 17β353β1201
Carmichael number:  19β43β409 19β199β271
Carmichael number:  23β199β353
Carmichael number:  29β113β1093 29β197β953
Carmichael number:  31β61β211 31β151β1171 31β181β331 31β271β601 31β991β15361
Carmichael number:  37β73β109 37β109β2017 37β613β1621
Carmichael number:  41β61β101 41β73β137 41β101β461 41β241β521 41β881β12041 41β1721β35281
Carmichael number:  43β127β211 43β211β337 43β271β5827 43β433β643 43β547β673 43β631β1597 43β3361β3907
Carmichael number:  47β1151β1933 47β3359β6073 47β3727β5153
Carmichael number:  53β79β599 53β157β521
Carmichael number:  59β1451β2089
Carmichael number:  61β181β1381 61β241β421 61β271β571 61β277β2113 61β421β12841 61β541β3001 61β661β2521 61β1301β19841 61β3361β4021
Carmichael number:  67β331β463 67β2311β51613
Carmichael number:  71β271β521 71β421β491 71β631β701 71β701β5531 71β911β9241
Carmichael number:  73β157β2293 73β379β523 73β601β21937 73β937β13681
Carmichael number:  79β2237β25247
Carmichael number:  83β2953β4019 83β6971β289297
Carmichael number:  89β353β617 89β617β3433 89β881β7129 89β4049β120121
Carmichael number:  97β193β1249 97β673β769 97β769β10657 97β1201β38833 97β2113β5857
Carmichael number:  101β151β251 101β1301β43801
Carmichael number:  103β647β1361 103β1123β6427 103β3571β183907 103β3877β4591 103β5407β185641
Carmichael number:  107β743β1061 107β3181β26183
Carmichael number:  109β163β379 109β229β4993 109β241β2389 109β379β919 109β433β2053 109β541β2269 109β1297β12853 109β1657β6229 109β1801β4789
Carmichael number:  113β337β449 113β827β18691 113β6833β85793 113β8737β22961
Carmichael number:  127β631β26713 127β659β1373 127β991β3313 127β2143β30241 127β16633β422479
Carmichael number:  131β491β4021 131β521β2731 131β571β1871 131β911β2341 131β1171β38351 131β1301β8971 131β4421β115831 131β5851β191621 131β17291β1132561
Carmichael number:  137β409β14009 137β953β5441 137β3673β20129 137β5441β11833
Carmichael number:  139β691β829 139β3359β66701 139β4003β92737 139β4969β8971 139β17389β21391
Carmichael number:  149β593β29453 149β1481β3109
Carmichael number:  151β211β541 151β331β3571 151β601β751 151β701β1451 151β751β8101 151β2551β192601 151β4951β53401
Carmichael number:  157β313β1093 157β937β6397 157β1093β15601 157β2017β28789 157β4993β11701 157β26053β409033
Carmichael number:  163β379β2377 163β487β811 163β811β1297 163β1297β42283 163β1621β37747
Carmichael number:  167β4483β34031 167β29383β490697
Carmichael number:  173β1291β31907 173β10321β255077 173β36809β155317
Carmichael number:  179β1069β4451 179β3739β7121 179β9257β57139 179β10859β485941
Carmichael number:  181β271β9811 181β541β3061 181β631β811 181β733β66337 181β1693β43777 181β2953β3637 181β3061β9721 181β5881β9421 181β11161β15661 181β16561β999181
Carmichael number:  191β421β431 191β571β15581 191β1901β7411 191β3877β56963 191β12541β342191
Carmichael number:  193β257β1601 193β401β11057 193β577β5569 193β1249β2593 193β2689β30529 193β7681β211777 193β61057β94273
Carmichael number:  199β397β4159 199β859β2311 199β937β20719 199β991β32869 199β4159β8713 199β8713β82567
Carmichael number:  211β281β491 211β421β631 211β631β66571 211β1051β9241 211β1741β4651 211β1831β4111 211β2311β54181 211β2521β3571 211β3221β35771 211β3571β9661 211β4019β11159 211β4201β98491 211β32341β70351 211β68041β127051
Carmichael number:  223β1777β23311 223β2221β5107 223β5107β37963 223β14653β79699 223β25087β1864801
Carmichael number:  227β1583β6781 227β2713β5651 227β7459β423299 227β21019β91757
Carmichael number:  229β457β2053 229β571β4219 229β1597β182857 229β2243β73379 229β7069β32377
Carmichael number:  233β5569β185369
Carmichael number:  239β409β3911 239β1429β1667 239β1667β6427 239β32369β234431
Carmichael number:  241β337β4513 241β433β52177 241β1201β2161 241β1361β4001 241β2161β3361 241β3361β32401 241β5521β110881 241β6481β780961 241β7321β588121 241β84961β181201
Carmichael number:  251β751β3251 251β2251β56501 251β3251β4001 251β4751β22501 251β16001β803251 251β22501β297251 251β31751β2656501
Carmichael number:  257β641β1153 257β769β49409 257β67073β3447553 257β78593β403969
Carmichael number:  263β787β2621 263β1049β3407 263β6551β12577 263β71527β1881161
Carmichael number:  269β1877β126229 269β4289β384581 269β10453β65393 269β15277β21977
Carmichael number:  271β541β811 271β811β2971 271β1171β7741 271β1801β16831 271β2161β65071 271β4591β4861 271β4861β77491 271β8191β1109881 271β8641β47791 271β9631β52201 271β10531β1426951 271β14797β1336663
Carmichael number:  277β829β7177 277β1381β1933 277β3313β9661 277β6073β31741 277β18493β88321 277β19597β36433
Carmichael number:  281β421β701 281β617β673 281β1009β4649 281β1321β23201 281β4201β9521 281β7121β26681 281β9521β13721 281β25621β84701
Carmichael number:  283β4231β598687 283β17203β58657
Carmichael number:  293β877β4673 293β1607β31391 293β3943β5987
Carmichael number:  307β613β919 307β919β141067 307β1531β3673 307β2143β13159 307β3673β225523 307β6427β246637 307β17443β153001 307β18973β1941571
Carmichael number:  311β1117β26723 311β1303β2357 311β2791β21701 311β3659β7069 311β23251β33791 311β26041β323951 311β28211β165541 311β44641β52391
Carmichael number:  313β521β23297 313β937β58657 313β1093β6709 313β1249β55849 313β3433β38377 313β3793β395737 313β5449β12097 313β6577β8761 313β7177β70201 313β9049β472057 313β12637β359581 313β49297β5143321 313β51481β947857 313β66457β184081 313β129169β400297
Carmichael number:  317β18013β41081 317β104281β2542853
Carmichael number:  331β661β991 331β947β24113 331β991β4621 331β1321β17491 331β2311β12541 331β2971β49171 331β3331β551281 331β4051β18121 331β4621β14851 331β37181β1758131 331β37951β897271 331β41141β316691
Carmichael number:  337β421β47293 337β449β21617 337β673β1009 337β953β8681 337β1009β3697 337β1597β12517 337β2473β11113 337β3697β12097 337β16369β1379089 337β19489β597073 337β35393β40433 337β58129β2176609
Carmichael number:  347β3461β92383 347β4153β29411
Carmichael number:  349β929β7541 349β1741β2089 349β3191β12239 349β4177β20533 349β122149β21315001
Carmichael number:  353β617β19801 353β1153β5153 353β1321β66617 353β13217β77761 353β130241β2704417
Carmichael number:  359β43319β3887881 359β46183β592133
Carmichael number:  367β733β1831 367β1831β9883 367β5003β42701 367β9151β419803 367β24889β51607 367β28183β574621
Carmichael number:  373β1117β1861 373β1613β150413 373β5581β1040857 373β16741β81097 373β139501β26016937
Carmichael number:  379β631β9199 379β757β4159 379β2269β24571 379β2539β21871 379β6427β202987 379β9829β17011 379β10639β268813
Carmichael number:  383β33617β40111 383β38201β860647 383β74873β3186263
Carmichael number:  389β3299β6791
Carmichael number:  397β1783β4951 397β2971β51283 397β4357β8317 397β30097β56629 397β55837β852589 397β79201β10480933 397β99793β370261
Carmichael number:  401β641β2161 401β1201β1601 401β2161β216641 401β2801β9601 401β9521β19681 401β9601β70001 401β15601β18401 401β18401β567601 401β161201β32320801
Carmichael number:  409β2857β6529 409β6121β96289 409β6529β22441 409β7039β575791 409β35089β683401 409β36721β114649
Carmichael number:  419β15467β47653 419β22573β78167 419β47653β539639
Carmichael number:  421β631β11551 421β701β2381 421β3851β85331 421β7561β289381 421β9661β15121 421β13441β209581 421β18481β39901 421β20231β54251 421β35533β7479697 421β42589β208489 421β89041β12495421
Carmichael number:  431β1721β29671 431β1979β142159 431β8171β55901 431β13331β168991
Carmichael number:  433β937β11593 433β1297β2161 433β2161β16417 433β2593β48817 433β2953β21673 433β3457β6481 433β3697β55201 433β6481β87697
Carmichael number:  439β3067β673207 439β3943β45553 439β9199β2019181 439β10513β17959 439β64679β7098521 439β96799β14164921
Carmichael number:  443β1327β4421 443β2029β4967 443β74257β143651 443β102103β2380613 443β167077β236471 443β251057β889747
Carmichael number:  449β2689β3137 449β50849β4566241 449β145601β325249 449β202049β45360001
Carmichael number:  457β3877β93253 457β5701β8893 457β7297β32377 457β15733β19381 457β21433β163249 457β28729β55633 457β71593β2337001 457β73721β1203233 457β114001β1211593
Carmichael number:  461β691β1151 461β1013β38917 461β1381β159161 461β3541β23321 461β5981β24841 461β26681β4099981
Carmichael number:  463β2927β15401 463β6007β39733 463β214831β49733377 463β218527β10117801
Carmichael number:  467β141199β474389
Carmichael number:  479β57839β219881
Carmichael number:  487β1459β8263 487β1531β2683 487β1621β1783 487β1783β108541 487β8263β9721 487β12637β32563 487β17011β2761453 487β26731β110323 487β51517β69499
Carmichael number:  491β1471β10781 491β6959β569479 491β16661β154351 491β41651β46061 491β122501β6683111 491β386611β637001
Carmichael number:  499β997β4483 499β10459β39841
Carmichael number:  503β5021β21587
Carmichael number:  509β3557β41149 509β7621β23369 509β11939β110491 509β86869β11054081
Carmichael number:  521β1301β8581 521β21841β41081
Carmichael number:  523β1567β163909 523β6091β1592797 523β9397β140419 523β15661β481807 523β38629β69427 523β155557β1114471 523β193663β462493
Carmichael number:  541β811β3511 541β1621β7561 541β6661β257401 541β7561β54541 541β12421β197641 541β16561β814501
Carmichael number:  547β1093β2731 547β2731β6553 547β6553β35491 547β7333β235951 547β26209β186187 547β52963β827737 547β158341β2624623
Carmichael number:  557β1669β42257 557β38921β7226333
Carmichael number:  563β28663β329333
Carmichael number:  569β2273β117577 569β13633β1108169 569β17609β25561 569β21017β37489 569β22153β787817
Carmichael number:  571β661β16411 571β2281β2851 571β2851β13681 571β6841β43891 571β13681β1562371 571β65323β18649717
Carmichael number:  577β757β39709 577β1153β6337 577β5569β100417 577β6337β26497 577β20161β646273 577β32833β37441 577β53857β181729 577β79777β86689 577β339841β15083713 577β559297β819073
Carmichael number:  587β5861β33403 587β9377β54499 587β12893β36919 587β49811β3654883
Carmichael number:  593β21017β31081 593β35521β3009137 593β176417β34871761
Carmichael number:  599β2393β84319 599β120199β17999801 599β179999β35939801 599β266111β547769 599β368369β12979591
Carmichael number:  601β1201β1801 601β1801β541201 601β3001β200401 601β3121β38281 601β3301β5101 601β4201β4801 601β4801β412201 601β5101β278701 601β6151β7951 601β9001β386401 601β19801β28201 601β52201β3921601 601β99901β923701
Carmichael number:  607β1213β9091 607β4243β1287751 607β21817β322999 607β24847β1885267 607β61813β7504099 607β186649β12588439 607β370873β45023983 607β373903β22695913
Carmichael number:  613β919β2143 613β1021β312937 613β1327β73951 613β1429β23053 613β2857β17341 613β7549β87313 613β9181β2813977 613β12241β111997 613β51817β213181 613β246637β783361 613β364753β386173
Carmichael number:  617β661β1013 617β8009β705937 617β16633β120737 617β29569β2606297 617β59753β81929 617β73613β133981 617β129361β6139673 617β137369β1629937 617β383153β47281081
Carmichael number:  619β1237β4327 619β2267β26987 619β5563β1721749 619β28429β703903 619β53149β56239 619β92083β452377 619β398611β9490009
Carmichael number:  631β1471β46411 631β5881β90511 631β26209β82279 631β32831β67481
Carmichael number:  641β4481β7681 641β12161β26881 641β17921β370561 641β19841β176641
Carmichael number:  643β107857β2391451
Carmichael number:  647β4523β19381 647β64601β75583 647β188633β532951 647β444449β7013623
Carmichael number:  653β13367β2909551 653β176041β732197
Carmichael number:  659β2633β5923 659β23689β624443 659β27919β34781 659β30269β92779 659β73039β6876101 659β92779β1329161
Carmichael number:  661β991β131011 661β1321β4621 661β2131β4231 661β3191β6491 661β3301β12541 661β4621β763621 661β5281β81181 661β22111β1623931 661β22441β95701 661β138821β152681
Carmichael number:  673β1009β14449 673β2017β3361 673β3361β12097 673β13441β1292257 673β40801β155137 673β231841β9178177
Carmichael number:  677β2029β85853 677β4733β1602121 677β6761β25013 677β45293β511057
Carmichael number:  683β8867β16369 683β11161β206027 683β15749β32303 683β42967β2934647 683β94117β9183131
Carmichael number:  691β7591β2622691 691β16561β2288731 691β31051β71761 691β34501β2648911 691β69691β3009781 691β743131β1330321
Carmichael number:  701β2801β10501 701β3701β1297201 701β3851β899851 701β6301β7001 701β18401β58901 701β41651β2245951 701β44101β170801 701β46901β319201 701β52501β296801 701β53201β632101
Carmichael number:  709β4957β12037 709β7789β16993 709β9677β21713 709β36109β5120257 709β210277β819157
Carmichael number:  719β97649β190271
Carmichael number:  727β1453β2179 727β2179β792067 727β2663β193601 727β3631β8713 727β4423β321553 727β176903β32152121 727β308551β1823713 727β651223β2784937
Carmichael number:  733β5857β84181 733β13177β47581 733β18301β789097 733β22571β2363507 733β25621β9390097 733β150427β1238911 733β271573β22118113 733β631717β3561913
Carmichael number:  739β821β4019 739β3691β454609 739β10333β2545363 739β62731β1783009 739β152029β1321759
Carmichael number:  743β6679β225569 743β6997β9011 743β596569β7266407
Carmichael number:  751β2251β10501 751β2851β237901 751β21751β181501 751β109751β649001 751β123001β1338751 751β153001β1767751 751β191251β10259251 751β318751β2418001
Carmichael number:  757β2017β18397 757β2269β858817 757β15121β3815533 757β27541β79273 757β32257β2219869 757β33013β59221 757β184843β633151 757β627481β6506893
Carmichael number:  761β2129β31769 761β2281β3041 761β3041β771401 761β6841β19001 761β8969β1137569 761β13681β101081 761β19001β1032841 761β41801β497041 761β230281β1184081 761β251941β339341 761β314641β497801
Carmichael number:  769β6529β9601 769β41729β697601
Carmichael number:  773β22003β122363 773β44777β47093
Carmichael number:  787β3931β9433 787β5503β45589 787β106373β3348623
Carmichael number:  797β2389β476009 797β3583β16319 797β5573β11941 797β21493β428249 797β58109β7718813 797β148853β859681
Carmichael number:  809β5657β9697 809β78781β176549 809β82013β22116173 809β176549β2197357 809β453289β1171601
Carmichael number:  811β1621β438211 811β4051β19441 811β4591β744661 811β6481β17011 811β19441β3153331 811β77761β1189891 811β86131β478441
Carmichael number:  821β1231β6971 821β15581β42641 821β137597β6275953
Carmichael number:  823β2467β4111 823β4111β23017 823β4933β9043 823β27127β637873 823β341953β31269703
Carmichael number:  827β2243β2833 827β4957β5783 827β24781β476603 827β101009β2880499 827β691363β57175721
Carmichael number:  829β1657β17389 829β9109β15733 829β10949β2269181 829β24841β1872109 829β140761β5556709
Carmichael number:  839β5867β223747
Carmichael number:  853β2557β4261 853β7669β594697 853β12781β5451097 853β17041β309277 853β19597β185737
Carmichael number:  857β6421β127973 857β10273β160073 857β95873β115561 857β796937β9229393
Carmichael number:  859β2861β3719 859β8581β9439 859β9439β150151 859β27457β66067 859β321751β1039039
Carmichael number:  863β24137β38791 863β28447β153437 863β38791β62927 863β56893β68099
Carmichael number:  877β1753β56941 877β3067β30223 877β6133β8761 877β24091β7042603 877β36793β6453493 877β263677β8894029
Carmichael number:  881β2861β840181 881β22441β57641 881β130241β16391761
Carmichael number:  883β2647β44101 883β8191β267877 883β11467β35281 883β15877β824671 883β16633β358219 883β21757β3842287 883β30871β134947 883β42337β216091 883β126127β161407 883β260191β114874327 883β403957β10808911 883β507151β531847
Carmichael number:  887β14177β50503
Carmichael number:  907β7853β16007 907β137713β24981139
Carmichael number:  911β2003β912367 911β9283β1208117 911β9311β55441 911β11831β898171 911β16381β28211 911β30941β4026751 911β55511β12642631 911β167441β204751 911β175631β2962961 911β185641β1551551 911β227501β2328691
Carmichael number:  919β8263β949213 919β15607β170749 919β60589β11136259 919β129439β569161 919β156979β321301 919β311203β2918323 919β877609β21797911
Carmichael number:  929β5569β23201 929β6961β35729 929β42689β1071841 929β139201β307169
Carmichael number:  937β1873β70201 937β6553β7489 937β7489β1002457 937β21529β3362113 937β38377β5993209 937β177841β820873
Carmichael number:  941β5171β23971 941β6581β8461 941β8461β361901 941β28201β102461 941β44651β4668511 941β209621β1133641 941β322891β701711 941β355321β1732421
Carmichael number:  947β29327β1983763 947β47129β299539 947β307451β10398433
Carmichael number:  953β2857β9521 953β5881β18257 953β17137β69497 953β52361β159937 953β159937β2771273
Carmichael number:  967β1289β25439 967β1933β4831 967β4831β11593 967β26083β5044453 967β62791β7589863 967β88873β1909783 967β156493β30265747
Carmichael number:  971β3881β753691 971β8731β44621 971β12611β3061321 971β110581β635351 971β142591β2387171 971β169751β648931 971β1324051β3263081
Carmichael number:  977β2441β794953 977β5857β12689 977β6833β39041 977β17569β41969 977β478241β155747153
Carmichael number:  983β3929β8839 983β8839β1241249 983β970217β190744663
Carmichael number:  991β4951β58411 991β10111β501001 991β16831β26731 991β56431β607861 991β99991β5215321 991β118801β206911 991β138403β336997 991β167311β312841 991β338581β890011 991β658351β1924561
Carmichael number:  997β1993β56773 997β8467β367027 997β12451β4137883 997β17929β130477 997β29383β450691 997β167329β15166093 997β1002973β99996409

ββββββββ  1038  Carmichael numbers found.


## Ruby

Works with: Ruby version 1.9
# Generate Charmichael Numbers require 'prime' Prime.each(61) do |p|  (2...p).each do |h3|    g = h3 + p    (1...g).each do |d|      next if (g*(p-1)) % d != 0 or (-p*p) % h3 != d % h3      q = 1 + ((p - 1) * g / d)      next unless q.prime?      r = 1 + (p * q / h3)      next unless r.prime? and (q * r) % (p - 1) == 1      puts "#{p} x #{q} x #{r}"     end  end  putsend
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


## Rust

 fn is_prime(n: i64) -> bool {    if n > 1 {        (2..((n / 2) + 1)).all(|x| n % x != 0)    } else {        false    }} // The modulo operator actually calculates the remainder.fn modulo(n: i64, m: i64) -> i64 {    ((n % m) + m) % m} fn carmichael(p1: i64) -> Vec<(i64, i64, i64)> {    let mut results = Vec::new();    if !is_prime(p1) {        return results;    }     for h3 in 2..p1 {        for d in 1..(h3 + p1) {            if (h3 + p1) * (p1 - 1) % d != 0 || modulo(-p1 * p1, h3) != d % h3 {                continue;            }             let p2 = 1 + ((p1 - 1) * (h3 + p1) / d);            if !is_prime(p2) {                continue;            }             let p3 = 1 + (p1 * p2 / h3);            if !is_prime(p3) || ((p2 * p3) % (p1 - 1) != 1) {                continue;            }             results.push((p1, p2, p3));        }    }     results} fn main() {    (1..62)        .filter(|&x| is_prime(x))        .map(carmichael)        .filter(|x| !x.is_empty())        .flat_map(|x| x)        .inspect(|x| println!("{:?}", x))        .count(); // Evaluate entire iterator}
Output:
(3, 11, 17)
(5, 29, 73)
(5, 17, 29)
(5, 13, 17)
.
.
.
(61, 661, 2521)
(61, 271, 571)
(61, 241, 421)
(61, 3361, 4021)


## Seed7

The function isPrime below is borrowed from the Seed7 algorithm collection.

$include "seed7_05.s7i"; const func boolean: isPrime (in integer: number) is func result var boolean: prime is FALSE; local var integer: upTo is 0; var integer: testNum is 3; begin if number = 2 then prime := TRUE; elsif odd(number) and number > 2 then upTo := sqrt(number); while number rem testNum <> 0 and testNum <= upTo do testNum +:= 2; end while; prime := testNum > upTo; end if; end func; const proc: main is func local var integer: p1 is 0; var integer: h3 is 0; var integer: g is 0; var integer: d is 0; var integer: p2 is 0; var integer: p3 is 0; begin for p1 range 2 to 61 do if isPrime(p1) then for h3 range 2 to p1 do g := h3 + p1; for d range 1 to pred(g) do if (g * pred(p1)) mod d = 0 and -p1 ** 2 mod h3 = d mod h3 then p2 := 1 + pred(p1) * g div d; if isPrime(p2) then p3 := 1 + p1 * p2 div h3; if isPrime(p3) and (p2 * p3) mod pred(p1) = 1 then writeln(p1 <& " * " <& p2 <& " * " <& p3 <& " = " <& p1*p2*p3); end if; end if; end if; end for; end for; end if; end for; end func; Output: 3 * 11 * 17 = 561 5 * 29 * 73 = 10585 5 * 17 * 29 = 2465 5 * 13 * 17 = 1105 7 * 19 * 67 = 8911 7 * 31 * 73 = 15841 7 * 13 * 31 = 2821 7 * 23 * 41 = 6601 7 * 73 * 103 = 52633 7 * 13 * 19 = 1729 13 * 61 * 397 = 314821 13 * 37 * 241 = 115921 13 * 97 * 421 = 530881 13 * 37 * 97 = 46657 13 * 37 * 61 = 29341 17 * 41 * 233 = 162401 17 * 353 * 1201 = 7207201 19 * 43 * 409 = 334153 19 * 199 * 271 = 1024651 23 * 199 * 353 = 1615681 29 * 113 * 1093 = 3581761 29 * 197 * 953 = 5444489 31 * 991 * 15361 = 471905281 31 * 61 * 631 = 1193221 31 * 151 * 1171 = 5481451 31 * 61 * 271 = 512461 31 * 61 * 211 = 399001 31 * 271 * 601 = 5049001 31 * 181 * 331 = 1857241 37 * 109 * 2017 = 8134561 37 * 73 * 541 = 1461241 37 * 613 * 1621 = 36765901 37 * 73 * 181 = 488881 37 * 73 * 109 = 294409 41 * 1721 * 35281 = 2489462641 41 * 881 * 12041 = 434932961 41 * 101 * 461 = 1909001 41 * 241 * 761 = 7519441 41 * 241 * 521 = 5148001 41 * 73 * 137 = 410041 41 * 61 * 101 = 252601 43 * 631 * 13567 = 368113411 43 * 271 * 5827 = 67902031 43 * 127 * 2731 = 14913991 43 * 127 * 1093 = 5968873 43 * 211 * 757 = 6868261 43 * 631 * 1597 = 43331401 43 * 127 * 211 = 1152271 43 * 211 * 337 = 3057601 43 * 433 * 643 = 11972017 43 * 547 * 673 = 15829633 43 * 3361 * 3907 = 564651361 47 * 3359 * 6073 = 958762729 47 * 1151 * 1933 = 104569501 47 * 3727 * 5153 = 902645857 53 * 157 * 2081 = 17316001 53 * 79 * 599 = 2508013 53 * 157 * 521 = 4335241 59 * 1451 * 2089 = 178837201 61 * 421 * 12841 = 329769721 61 * 181 * 5521 = 60957361 61 * 1301 * 19841 = 1574601601 61 * 277 * 2113 = 35703361 61 * 181 * 1381 = 15247621 61 * 541 * 3001 = 99036001 61 * 661 * 2521 = 101649241 61 * 271 * 571 = 9439201 61 * 241 * 421 = 6189121 61 * 3361 * 4021 = 824389441  ## Sidef Translation of: Perl func forprimes(a, b, callback) { for (a = (a-1 -> next_prime); a <= b; a.next_prime!) { callback(a) }} forprimes(3, 61, func(p) { for h3 in (2 .. p-1) { var ph3 = (p + h3) for d in (1 .. ph3-1) { ((-p * p) % h3) != (d % h3) && next ((p-1) * ph3) % d && next var q = 1+((p-1) * ph3 / d) q.is_prime || next var r = 1+((p*q - 1)/h3) r.is_prime || next (q*r) % (p-1) == 1 || next printf("%2d x %5d x %5d = %s\n",p,q,r, p*q*r) } }}) Output:  3 x 11 x 17 = 561 5 x 29 x 73 = 10585 5 x 17 x 29 = 2465 5 x 13 x 17 = 1105 ... full output is 69 lines ... 61 x 661 x 2521 = 101649241 61 x 271 x 571 = 9439201 61 x 241 x 421 = 6189121 61 x 3361 x 4021 = 824389441  ## Tcl Using the primality tester from the Miller-Rabin task... proc carmichael {limit {rounds 10}} { set carmichaels {} for {set p1 2} {$p1 <= $limit} {incr p1} { if {![miller_rabin$p1 $rounds]} continue for {set h3 2} {$h3 < $p1} {incr h3} { set g [expr {$h3 + $p1}] for {set d 1} {$d < $h3+$p1} {incr d} {		if {(($h3+$p1)*($p1-1))%$d != 0} continue		if {(-($p1**2))%$h3 != $d%$h3} continue 		set p2 [expr {1 + ($p1-1)*$g/$d}] if {![miller_rabin$p2 $rounds]} continue set p3 [expr {1 +$p1*$p2/$h3}]		if {![miller_rabin $p3$rounds]} continue 		if {($p2*$p3)%($p1-1) != 1} continue lappend carmichaels$p1 $p2$p3 [expr {$p1*$p2*$p3}] } } } return$carmichaels}

Demonstrating:

set results [carmichael 61 2]puts "[expr {[llength $results]/4}] Carmichael numbers found"foreach {p1 p2 p3 c}$results {    puts "$p1 x$p2 x $p3 =$c"}
Output:
69 Carmichael numbers found
3 x 11 x 17 = 561
5 x 29 x 73 = 10585
5 x 17 x 29 = 2465
5 x 13 x 17 = 1105
7 x 19 x 67 = 8911
7 x 31 x 73 = 15841
7 x 13 x 31 = 2821
7 x 23 x 41 = 6601
7 x 73 x 103 = 52633
7 x 13 x 19 = 1729
13 x 61 x 397 = 314821
13 x 37 x 241 = 115921
13 x 97 x 421 = 530881
13 x 37 x 97 = 46657
13 x 37 x 61 = 29341
17 x 41 x 233 = 162401
17 x 353 x 1201 = 7207201
19 x 43 x 409 = 334153
19 x 199 x 271 = 1024651
23 x 199 x 353 = 1615681
29 x 113 x 1093 = 3581761
29 x 197 x 953 = 5444489
31 x 991 x 15361 = 471905281
31 x 61 x 631 = 1193221
31 x 151 x 1171 = 5481451
31 x 61 x 271 = 512461
31 x 61 x 211 = 399001
31 x 271 x 601 = 5049001
31 x 181 x 331 = 1857241
37 x 109 x 2017 = 8134561
37 x 73 x 541 = 1461241
37 x 613 x 1621 = 36765901
37 x 73 x 181 = 488881
37 x 73 x 109 = 294409
41 x 1721 x 35281 = 2489462641
41 x 881 x 12041 = 434932961
41 x 101 x 461 = 1909001
41 x 241 x 761 = 7519441
41 x 241 x 521 = 5148001
41 x 73 x 137 = 410041
41 x 61 x 101 = 252601
43 x 631 x 13567 = 368113411
43 x 271 x 5827 = 67902031
43 x 127 x 2731 = 14913991
43 x 127 x 1093 = 5968873
43 x 211 x 757 = 6868261
43 x 631 x 1597 = 43331401
43 x 127 x 211 = 1152271
43 x 211 x 337 = 3057601
43 x 433 x 643 = 11972017
43 x 547 x 673 = 15829633
43 x 3361 x 3907 = 564651361
47 x 3359 x 6073 = 958762729
47 x 1151 x 1933 = 104569501
47 x 3727 x 5153 = 902645857
53 x 157 x 2081 = 17316001
53 x 79 x 599 = 2508013
53 x 157 x 521 = 4335241
59 x 1451 x 2089 = 178837201
61 x 421 x 12841 = 329769721
61 x 181 x 5521 = 60957361
61 x 1301 x 19841 = 1574601601
61 x 277 x 2113 = 35703361
61 x 181 x 1381 = 15247621
61 x 541 x 3001 = 99036001
61 x 661 x 2521 = 101649241
61 x 271 x 571 = 9439201
61 x 241 x 421 = 6189121
61 x 3361 x 4021 = 824389441


## zkl

Using the Miller-Rabin primality test in lib GMP.

var BN=Import("zklBigNum"), bi=BN(0); // gonna recycle biprimes:=T(2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,59,61);var p2,p3;cs:=[[(p1,h3,d); primes; { [2..p1 - 1] }; // list comprehension      { [1..h3 + p1 - 1] },	{ ((h3 + p1)*(p1 - 1)%d == 0 and ((-p1*p1):mod(_,h3) == d%h3)) },//guard	{ (p2=1 + (p1 - 1)*(h3 + p1)/d):bi.set(_).probablyPrime() },//guard	{ (p3=1 + (p1*p2/h3)):bi.set(_).probablyPrime() },	 //guard	{ 1==(p2*p3)%(p1 - 1) };				 //guard   { T(p1,p2,p3) }  // return list of three primes in Carmichael number]];fcn mod(a,b) { m:=a%b; if(m<0) m+b else m }
cs.len().println(" Carmichael numbers found:");cs.pump(Console.println,fcn([(p1,p2,p3)]){   "%2d * %4d * %5d = %d".fmt(p1,p2,p3,p1*p2*p3) });
Output:
69 Carmichael numbers found:
3 *   11 *    17 = 561
5 *   29 *    73 = 10585
5 *   17 *    29 = 2465
5 *   13 *    17 = 1105
7 *   19 *    67 = 8911
...
61 *  181 *  1381 = 15247621
61 *  541 *  3001 = 99036001
61 *  661 *  2521 = 101649241
61 *  271 *   571 = 9439201
61 *  241 *   421 = 6189121
61 * 3361 *  4021 = 824389441


## ZX Spectrum Basic

Translation of: C
10 FOR p=2 TO 6120 LET n=p: GO SUB 100030 IF NOT n THEN GO TO 20040 FOR h=1 TO p-150 FOR d=1 TO h-1+p60 IF NOT (FN m((h+p)*(p-1),d)=0 AND FN w(-p*p,h)=FN m(d,h)) THEN GO TO 18070 LET q=INT (1+((p-1)*(h+p)/d))80 LET n=q: GO SUB 100090 IF NOT n THEN GO TO 180100 LET r=INT (1+(p*q/h))110 LET n=r: GO SUB 1000120 IF (NOT n) OR ((FN m((q*r),(p-1))<>1)) THEN GO TO 180130 PRINT p;" ";q;" ";r180 NEXT d190 NEXT h200 NEXT p210 STOP 1000 IF n<4 THEN LET n=(n>1): RETURN 1010 IF (NOT FN m(n,2)) OR (NOT FN m(n,3)) THEN LET n=0: RETURN 1020 LET i=51030 IF NOT ((i*i)<=n) THEN LET n=1: RETURN 1040 IF (NOT FN m(n,i)) OR NOT FN m(n,(i+2)) THEN LET n=0: RETURN 1050 LET i=i+61060 GO TO 10302000 DEF FN m(a,b)=a-(INT (a/b)*b): REM Mod function2010 DEF FN w(a,b)=FN m(FN m(a,b)+b,b): REM Mod function modified