Fermat pseudoprimes

A Fermat pseudoprime is a positive composite integer that passes the Fermat primality test.

Fermat pseudoprimes is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Fermat's little theorem states that if p is prime and a is coprime to p, then ap−1 − 1 is divisible by p.

For an integer a > 1, if a composite integer x evenly divides ax−1 − 1, then x is called a Fermat pseudoprime to base a.

Fermat pseudoprimes to base 2 are sometimes called Sarrus numbers or Poulet numbers.

Fermat pseudoprimes can be found to any positive integer base. When using a base integer a = 1, this method returns all composite numbers.


Task

For base integers a of 1 through 20:

  • Find the count of pseudoprimes up to and including 12,000.
  • Show the first 20 pseudoprimes.


Stretch
  • Extend the count threshold out to 25,000, 50,000 or higher.


See also


Related


ALGOL 68Edit

BEGIN # find some Fermat pseudo primes: x is a Fermat pseudoprime in base a   #
      # if a^(x-1) - 1 is divisible by x and x is not prime                   #

    # iterative Greatest Common Divisor routine, returns the gcd of m and n   #
    PROC gcd = ( INT m, n )INT:
         BEGIN
            INT a := ABS m, b := ABS n;
            WHILE b /= 0 DO
                INT new a = b;
                b        := a MOD b;
                a        := new a
            OD;
            a
        END # gcd # ;
    [ 0 : 50 000 ]BOOL prime;
    prime[ 0 ] := prime[ 1 ] := FALSE;
    prime[ 2 ] := TRUE;
    FOR i FROM 3 BY 2 TO UPB prime DO prime[ i ] := TRUE  OD;
    FOR i FROM 4 BY 2 TO UPB prime DO prime[ i ] := FALSE OD;
    FOR i FROM 3 BY 2 TO ENTIER sqrt( UPB prime ) DO
        IF prime[ i ] THEN
            FOR s FROM i * i BY i + i TO UPB prime DO prime[ s ] := FALSE OD
        FI
    OD;
    FOR base FROM 1 TO 20 DO
        [ 1 : 20 ]INT fp;
        INT count := 0;
        # x from 3 as 1 is neither composite nor prime and 2 is prime         #
        IF base = 1 THEN
            # 1^(x-1) is 1 for all x, so all composites are                   #
            # fermat pseudo primes in base 1                                  #
            FOR x FROM 3 TO UPB prime DO
                IF NOT prime[ x ] THEN
                    IF ( count +:= 1 ) <= UPB fp THEN
                        fp[ count ] := x
                    FI
                FI
            OD
        ELSE
            # must test base^(x-1) mod x                                      #
            FOR x FROM 3 BY IF ODD base THEN 1 ELSE 2 FI TO UPB prime DO
                IF NOT prime[ x ] THEN
                    IF gcd( x, base ) = 1 THEN
                        # have a composite x co-prime to the base             #
                        INT base to x minus 1 mod x := 1;
                        FOR p TO x - 1
                        WHILE
                            base to x minus 1 mod x *:= base MODAB x;
                            IF base to x minus 1 mod x = 1
                            THEN # the power is now 1 mod x,                  #
                                 # if the power divides x - 1 then base^(x-1) #
                                 # will be 1 mod x otherwise it won't be      #
                                 IF ( x - 1 ) MOD p /= 0 THEN
                                     base to x minus 1 mod x := 0
                                 FI;
                                 FALSE
                            ELSE TRUE # keep checking the powers              #
                            FI
                        DO SKIP OD;
                        IF base to x minus 1 mod x = 1 THEN
                            # have a composite x that divides base^(x-1)-1    #
                            IF ( count +:= 1 ) <= UPB fp THEN
                                fp[ count ] := x
                            FI
                        FI
                    FI
                FI
            OD
        FI;
        print( ( "base ", whole( base, -2 ), " to ", whole( UPB prime, 0 )
               , " total: ", whole( count, -6 ), ", first: "
               )
             );
        FOR f FROM LWB fp TO IF count < UPB fp THEN count ELSE UPB fp FI DO
            print( ( " ", whole( fp[ f ], 0 ) ) )
        OD;
        print( ( newline ) )
    OD
END
Output:
base  1 to 50000 total:  44866, first:  4 6 8 9 10 12 14 15 16 18 20 21 22 24 25 26 27 28 30 32
base  2 to 50000 total:     55, first:  341 561 645 1105 1387 1729 1905 2047 2465 2701 2821 3277 4033 4369 4371 4681 5461 6601 7957 8321
base  3 to 50000 total:     53, first:  91 121 286 671 703 949 1105 1541 1729 1891 2465 2665 2701 2821 3281 3367 3751 4961 5551 6601
base  4 to 50000 total:    111, first:  15 85 91 341 435 451 561 645 703 1105 1247 1271 1387 1581 1695 1729 1891 1905 2047 2071
base  5 to 50000 total:     54, first:  4 124 217 561 781 1541 1729 1891 2821 4123 5461 5611 5662 5731 6601 7449 7813 8029 8911 9881
base  6 to 50000 total:     74, first:  35 185 217 301 481 1105 1111 1261 1333 1729 2465 2701 2821 3421 3565 3589 3913 4123 4495 5713
base  7 to 50000 total:     49, first:  6 25 325 561 703 817 1105 1825 2101 2353 2465 3277 4525 4825 6697 8321 10225 10585 10621 11041
base  8 to 50000 total:    150, first:  9 21 45 63 65 105 117 133 153 231 273 341 481 511 561 585 645 651 861 949
base  9 to 50000 total:    113, first:  4 8 28 52 91 121 205 286 364 511 532 616 671 697 703 946 949 1036 1105 1288
base 10 to 50000 total:     65, first:  9 33 91 99 259 451 481 561 657 703 909 1233 1729 2409 2821 2981 3333 3367 4141 4187
base 11 to 50000 total:     61, first:  10 15 70 133 190 259 305 481 645 703 793 1105 1330 1729 2047 2257 2465 2821 4577 4921
base 12 to 50000 total:     91, first:  65 91 133 143 145 247 377 385 703 1045 1099 1105 1649 1729 1885 1891 2041 2233 2465 2701
base 13 to 50000 total:     68, first:  4 6 12 21 85 105 231 244 276 357 427 561 1099 1785 1891 2465 2806 3605 5028 5149
base 14 to 50000 total:     69, first:  15 39 65 195 481 561 781 793 841 985 1105 1111 1541 1891 2257 2465 2561 2665 2743 3277
base 15 to 50000 total:     42, first:  14 341 742 946 1477 1541 1687 1729 1891 1921 2821 3133 3277 4187 6541 6601 7471 8701 8911 9073
base 16 to 50000 total:    145, first:  15 51 85 91 255 341 435 451 561 595 645 703 1105 1247 1261 1271 1285 1387 1581 1687
base 17 to 50000 total:     63, first:  4 8 9 16 45 91 145 261 781 1111 1228 1305 1729 1885 2149 2821 3991 4005 4033 4187
base 18 to 50000 total:     98, first:  25 49 65 85 133 221 323 325 343 425 451 637 931 1105 1225 1369 1387 1649 1729 1921
base 19 to 50000 total:     93, first:  6 9 15 18 45 49 153 169 343 561 637 889 905 906 1035 1105 1629 1661 1849 1891
base 20 to 50000 total:     66, first:  21 57 133 231 399 561 671 861 889 1281 1653 1729 1891 2059 2413 2501 2761 2821 2947 3059

C++Edit

Translation of: Wren
#include <cstdint>
#include <iomanip>
#include <iostream>

uint64_t modpow(uint64_t base, uint64_t exp, uint64_t mod) {
    if (mod == 1)
        return 0;
    uint64_t result = 1;
    base %= mod;
    for (; exp > 0; exp >>= 1) {
        if ((exp & 1) == 1)
            result = (result * base) % mod;
        base = (base * base) % mod;
    }
    return result;
}

bool is_prime(uint64_t n) {
    if (n < 2)
        return false;
    if (n % 2 == 0)
        return n == 2;
    if (n % 3 == 0)
        return n == 3;
    for (uint64_t p = 5; p * p <= n; p += 4) {
        if (n % p == 0)
            return false;
        p += 2;
        if (n % p == 0)
            return false;
    }
    return true;
}

bool is_fermat_pseudoprime(uint64_t a, uint64_t x) {
    return !is_prime(x) && modpow(a, x - 1, x) == 1;
}

int main() {
    std::cout << "First 20 Fermat pseudoprimes:\n";
    for (uint64_t a = 1; a <= 20; ++a) {
        std::cout << "Base " << std::setw(2) << a << ": ";
        int count = 0;
        for (uint64_t x = 4; count < 20; ++x) {
            if (is_fermat_pseudoprime(a, x)) {
                ++count;
                std::cout << std::setw(5) << x << ' ';
            }
        }
        std::cout << '\n';
    }

    const uint64_t limits[] = {12000, 25000, 50000, 100000};
    std::cout << "\nCount <= ";
    for (uint64_t limit : limits) {
        std::cout << std::setw(6) << limit << ' ';
    }
    std::cout << "\n------------------------------------\n";
    for (uint64_t a = 1; a <= 20; ++a) {
        std::cout << "Base " << std::setw(2) << a << ": ";
        int count = 0;
        uint64_t x = 4;
        for (uint64_t limit : limits) {
            for (; x <= limit; ++x) {
                if (is_fermat_pseudoprime(a, x))
                    ++count;
            }
            std::cout << std::setw(6) << count << ' ';
        }
        std::cout << '\n';
    }
}
Output:
First 20 Fermat pseudoprimes:
Base  1:     4     6     8     9    10    12    14    15    16    18    20    21    22    24    25    26    27    28    30    32 
Base  2:   341   561   645  1105  1387  1729  1905  2047  2465  2701  2821  3277  4033  4369  4371  4681  5461  6601  7957  8321 
Base  3:    91   121   286   671   703   949  1105  1541  1729  1891  2465  2665  2701  2821  3281  3367  3751  4961  5551  6601 
Base  4:    15    85    91   341   435   451   561   645   703  1105  1247  1271  1387  1581  1695  1729  1891  1905  2047  2071 
Base  5:     4   124   217   561   781  1541  1729  1891  2821  4123  5461  5611  5662  5731  6601  7449  7813  8029  8911  9881 
Base  6:    35   185   217   301   481  1105  1111  1261  1333  1729  2465  2701  2821  3421  3565  3589  3913  4123  4495  5713 
Base  7:     6    25   325   561   703   817  1105  1825  2101  2353  2465  3277  4525  4825  6697  8321 10225 10585 10621 11041 
Base  8:     9    21    45    63    65   105   117   133   153   231   273   341   481   511   561   585   645   651   861   949 
Base  9:     4     8    28    52    91   121   205   286   364   511   532   616   671   697   703   946   949  1036  1105  1288 
Base 10:     9    33    91    99   259   451   481   561   657   703   909  1233  1729  2409  2821  2981  3333  3367  4141  4187 
Base 11:    10    15    70   133   190   259   305   481   645   703   793  1105  1330  1729  2047  2257  2465  2821  4577  4921 
Base 12:    65    91   133   143   145   247   377   385   703  1045  1099  1105  1649  1729  1885  1891  2041  2233  2465  2701 
Base 13:     4     6    12    21    85   105   231   244   276   357   427   561  1099  1785  1891  2465  2806  3605  5028  5149 
Base 14:    15    39    65   195   481   561   781   793   841   985  1105  1111  1541  1891  2257  2465  2561  2665  2743  3277 
Base 15:    14   341   742   946  1477  1541  1687  1729  1891  1921  2821  3133  3277  4187  6541  6601  7471  8701  8911  9073 
Base 16:    15    51    85    91   255   341   435   451   561   595   645   703  1105  1247  1261  1271  1285  1387  1581  1687 
Base 17:     4     8     9    16    45    91   145   261   781  1111  1228  1305  1729  1885  2149  2821  3991  4005  4033  4187 
Base 18:    25    49    65    85   133   221   323   325   343   425   451   637   931  1105  1225  1369  1387  1649  1729  1921 
Base 19:     6     9    15    18    45    49   153   169   343   561   637   889   905   906  1035  1105  1629  1661  1849  1891 
Base 20:    21    57   133   231   399   561   671   861   889  1281  1653  1729  1891  2059  2413  2501  2761  2821  2947  3059 

Count <=  12000  25000  50000 100000 
------------------------------------
Base  1:  10561  22237  44866  90407 
Base  2:     25     38     55     78 
Base  3:     25     38     53     78 
Base  4:     50     75    111    153 
Base  5:     22     35     54     73 
Base  6:     31     46     74    104 
Base  7:     21     32     49     73 
Base  8:     76    110    150    218 
Base  9:     55     83    113    164 
Base 10:     35     53     65     90 
Base 11:     30     41     61     89 
Base 12:     35     60     91    127 
Base 13:     31     51     68     91 
Base 14:     33     51     69     96 
Base 15:     22     31     42     70 
Base 16:     69    102    145    200 
Base 17:     31     44     63     83 
Base 18:     46     69     98    134 
Base 19:     48     70     93    121 
Base 20:     35     49     66    101 

F#Edit

This task uses Extensible Prime Generator (F#)

// Fermat pseudoprimes. Nigel Galloway: August 17th., 2022
let fp(a:int)=let a=bigint a in primesI()|>Seq.pairwise|>Seq.collect(fun(n,g)->seq{for n in n+1I..g-1I do if bigint.ModPow(a,n-1I,n)=1I then yield n})
{1..20}|>Seq.iter(fun n->printf $"Base %2d{n} - Up to 50000: %5d{fp n|>Seq.takeWhile((>=)50000I)|>Seq.length} First 20: ("; fp n|>Seq.take 20|>Seq.iter(printf "%A "); printfn ")")
Output:
Base  1 - Up to 50000: 44866 First 20: (4 6 8 9 10 12 14 15 16 18 20 21 22 24 25 26 27 28 30 32 )
Base  2 - Up to 50000:    55 First 20: (341 561 645 1105 1387 1729 1905 2047 2465 2701 2821 3277 4033 4369 4371 4681 5461 6601 7957 8321 )
Base  3 - Up to 50000:    53 First 20: (91 121 286 671 703 949 1105 1541 1729 1891 2465 2665 2701 2821 3281 3367 3751 4961 5551 6601 )
Base  4 - Up to 50000:   111 First 20: (15 85 91 341 435 451 561 645 703 1105 1247 1271 1387 1581 1695 1729 1891 1905 2047 2071 )
Base  5 - Up to 50000:    54 First 20: (4 124 217 561 781 1541 1729 1891 2821 4123 5461 5611 5662 5731 6601 7449 7813 8029 8911 9881 )
Base  6 - Up to 50000:    74 First 20: (35 185 217 301 481 1105 1111 1261 1333 1729 2465 2701 2821 3421 3565 3589 3913 4123 4495 5713 )
Base  7 - Up to 50000:    49 First 20: (6 25 325 561 703 817 1105 1825 2101 2353 2465 3277 4525 4825 6697 8321 10225 10585 10621 11041 )
Base  8 - Up to 50000:   150 First 20: (9 21 45 63 65 105 117 133 153 231 273 341 481 511 561 585 645 651 861 949 )
Base  9 - Up to 50000:   113 First 20: (4 8 28 52 91 121 205 286 364 511 532 616 671 697 703 946 949 1036 1105 1288 )
Base 10 - Up to 50000:    65 First 20: (9 33 91 99 259 451 481 561 657 703 909 1233 1729 2409 2821 2981 3333 3367 4141 4187 )
Base 11 - Up to 50000:    61 First 20: (10 15 70 133 190 259 305 481 645 703 793 1105 1330 1729 2047 2257 2465 2821 4577 4921 )
Base 12 - Up to 50000:    91 First 20: (65 91 133 143 145 247 377 385 703 1045 1099 1105 1649 1729 1885 1891 2041 2233 2465 2701 )
Base 13 - Up to 50000:    68 First 20: (4 6 12 21 85 105 231 244 276 357 427 561 1099 1785 1891 2465 2806 3605 5028 5149 )
Base 14 - Up to 50000:    69 First 20: (15 39 65 195 481 561 781 793 841 985 1105 1111 1541 1891 2257 2465 2561 2665 2743 3277 )
Base 15 - Up to 50000:    42 First 20: (14 341 742 946 1477 1541 1687 1729 1891 1921 2821 3133 3277 4187 6541 6601 7471 8701 8911 9073 )
Base 16 - Up to 50000:   145 First 20: (15 51 85 91 255 341 435 451 561 595 645 703 1105 1247 1261 1271 1285 1387 1581 1687 )
Base 17 - Up to 50000:    63 First 20: (4 8 9 16 45 91 145 261 781 1111 1228 1305 1729 1885 2149 2821 3991 4005 4033 4187 )
Base 18 - Up to 50000:    98 First 20: (25 49 65 85 133 221 323 325 343 425 451 637 931 1105 1225 1369 1387 1649 1729 1921 )
Base 19 - Up to 50000:    93 First 20: (6 9 15 18 45 49 153 169 343 561 637 889 905 906 1035 1105 1629 1661 1849 1891 )
Base 20 - Up to 50000:    66 First 20: (21 57 133 231 399 561 671 861 889 1281 1653 1729 1891 2059 2413 2501 2761 2821 2947 3059 )
Real: 00:00:00.632

JuliaEdit

using Primes

ispseudo(n, base) = !isprime(n) && BigInt(base)^(n - 1) % n == 1

for b in 1:20
    pseudos = filter(n -> ispseudo(n, b), 1:50000)
    println("Base ", lpad(b, 2), " up to 50000: ", lpad(length(pseudos), 5), "  First 20: ", pseudos[1:20])
end
Output:
Base  1 up to 50000: 44866  First 20: [4, 6, 8, 9, 10, 12, 14, 15, 16, 18, 20, 21, 22, 24, 25, 26, 27, 28, 30, 32]
Base  2 up to 50000:    55  First 20: [341, 561, 645, 1105, 1387, 1729, 1905, 2047, 2465, 2701, 2821, 3277, 4033, 4369, 4371, 4681, 5461, 6601, 7957, 8321]      
Base  3 up to 50000:    53  First 20: [91, 121, 286, 671, 703, 949, 1105, 1541, 1729, 1891, 2465, 2665, 2701, 2821, 3281, 3367, 3751, 4961, 5551, 6601]
Base  4 up to 50000:   111  First 20: [15, 85, 91, 341, 435, 451, 561, 645, 703, 1105, 1247, 1271, 1387, 1581, 1695, 1729, 1891, 1905, 2047, 2071]
Base  5 up to 50000:    54  First 20: [4, 124, 217, 561, 781, 1541, 1729, 1891, 2821, 4123, 5461, 5611, 5662, 5731, 6601, 7449, 7813, 8029, 8911, 9881]
Base  6 up to 50000:    74  First 20: [35, 185, 217, 301, 481, 1105, 1111, 1261, 1333, 1729, 2465, 2701, 2821, 3421, 3565, 3589, 3913, 4123, 4495, 5713]
Base  7 up to 50000:    49  First 20: [6, 25, 325, 561, 703, 817, 1105, 1825, 2101, 2353, 2465, 3277, 4525, 4825, 6697, 8321, 10225, 10585, 10621, 11041]
Base  8 up to 50000:   150  First 20: [9, 21, 45, 63, 65, 105, 117, 133, 153, 231, 273, 341, 481, 511, 561, 585, 645, 651, 861, 949]
Base  9 up to 50000:   113  First 20: [4, 8, 28, 52, 91, 121, 205, 286, 364, 511, 532, 616, 671, 697, 703, 946, 949, 1036, 1105, 1288]
Base 10 up to 50000:    65  First 20: [9, 33, 91, 99, 259, 451, 481, 561, 657, 703, 909, 1233, 1729, 2409, 2821, 2981, 3333, 3367, 4141, 4187]
Base 11 up to 50000:    61  First 20: [10, 15, 70, 133, 190, 259, 305, 481, 645, 703, 793, 1105, 1330, 1729, 2047, 2257, 2465, 2821, 4577, 4921]
Base 12 up to 50000:    91  First 20: [65, 91, 133, 143, 145, 247, 377, 385, 703, 1045, 1099, 1105, 1649, 1729, 1885, 1891, 2041, 2233, 2465, 2701]
Base 13 up to 50000:    68  First 20: [4, 6, 12, 21, 85, 105, 231, 244, 276, 357, 427, 561, 1099, 1785, 1891, 2465, 2806, 3605, 5028, 5149]
Base 14 up to 50000:    69  First 20: [15, 39, 65, 195, 481, 561, 781, 793, 841, 985, 1105, 1111, 1541, 1891, 2257, 2465, 2561, 2665, 2743, 3277]
Base 15 up to 50000:    42  First 20: [14, 341, 742, 946, 1477, 1541, 1687, 1729, 1891, 1921, 2821, 3133, 3277, 4187, 6541, 6601, 7471, 8701, 8911, 9073]
Base 16 up to 50000:   145  First 20: [15, 51, 85, 91, 255, 341, 435, 451, 561, 595, 645, 703, 1105, 1247, 1261, 1271, 1285, 1387, 1581, 1687]
Base 17 up to 50000:    63  First 20: [4, 8, 9, 16, 45, 91, 145, 261, 781, 1111, 1228, 1305, 1729, 1885, 2149, 2821, 3991, 4005, 4033, 4187]
Base 18 up to 50000:    98  First 20: [25, 49, 65, 85, 133, 221, 323, 325, 343, 425, 451, 637, 931, 1105, 1225, 1369, 1387, 1649, 1729, 1921]
Base 19 up to 50000:    93  First 20: [6, 9, 15, 18, 45, 49, 153, 169, 343, 561, 637, 889, 905, 906, 1035, 1105, 1629, 1661, 1849, 1891]
Base 20 up to 50000:    66 First 20: [21, 57, 133, 231, 399, 561, 671, 861, 889, 1281, 1653, 1729, 1891, 2059, 2413, 2501, 2761, 2821, 2947, 3059]

PerlEdit

Library: ntheory
use v5.36;
use ntheory 'is_prime';

sub expmod ($a, $b, $n) {
    my $c = 1;
    do {
        ($c *= $a) %= $n if $b % 2;
        ($a *= $a) %= $n;
    } while ($b = int $b/2);
    $c
}

my $threshold = 50000;
say "For each base: # of terms up to $threshold, first 20 displayed";
for my $b (1..20) {
    my @pseudo = grep { !is_prime($_) && (1 == expmod $b, $_ - 1, $_) } 1..$threshold;
    printf "base %2d: %5d (%s)\n", $b, $#pseudo, join ' ', @pseudo[1..20];
}
Output:
base  1: 44866 (4 6 8 9 10 12 14 15 16 18 20 21 22 24 25 26 27 28 30 32)
base  2:    55 (341 561 645 1105 1387 1729 1905 2047 2465 2701 2821 3277 4033 4369 4371 4681 5461 6601 7957 8321)
base  3:    53 (91 121 286 671 703 949 1105 1541 1729 1891 2465 2665 2701 2821 3281 3367 3751 4961 5551 6601)
base  4:   111 (15 85 91 341 435 451 561 645 703 1105 1247 1271 1387 1581 1695 1729 1891 1905 2047 2071)
base  5:    54 (4 124 217 561 781 1541 1729 1891 2821 4123 5461 5611 5662 5731 6601 7449 7813 8029 8911 9881)
base  6:    74 (35 185 217 301 481 1105 1111 1261 1333 1729 2465 2701 2821 3421 3565 3589 3913 4123 4495 5713)
base  7:    49 (6 25 325 561 703 817 1105 1825 2101 2353 2465 3277 4525 4825 6697 8321 10225 10585 10621 11041)
base  8:   150 (9 21 45 63 65 105 117 133 153 231 273 341 481 511 561 585 645 651 861 949)
base  9:   113 (4 8 28 52 91 121 205 286 364 511 532 616 671 697 703 946 949 1036 1105 1288)
base 10:    65 (9 33 91 99 259 451 481 561 657 703 909 1233 1729 2409 2821 2981 3333 3367 4141 4187)
base 11:    61 (10 15 70 133 190 259 305 481 645 703 793 1105 1330 1729 2047 2257 2465 2821 4577 4921)
base 12:    91 (65 91 133 143 145 247 377 385 703 1045 1099 1105 1649 1729 1885 1891 2041 2233 2465 2701)
base 13:    68 (4 6 12 21 85 105 231 244 276 357 427 561 1099 1785 1891 2465 2806 3605 5028 5149)
base 14:    69 (15 39 65 195 481 561 781 793 841 985 1105 1111 1541 1891 2257 2465 2561 2665 2743 3277)
base 15:    42 (14 341 742 946 1477 1541 1687 1729 1891 1921 2821 3133 3277 4187 6541 6601 7471 8701 8911 9073)
base 16:   145 (15 51 85 91 255 341 435 451 561 595 645 703 1105 1247 1261 1271 1285 1387 1581 1687)
base 17:    63 (4 8 9 16 45 91 145 261 781 1111 1228 1305 1729 1885 2149 2821 3991 4005 4033 4187)
base 18:    98 (25 49 65 85 133 221 323 325 343 425 451 637 931 1105 1225 1369 1387 1649 1729 1921)
base 19:    93 (6 9 15 18 45 49 153 169 343 561 637 889 905 906 1035 1105 1629 1661 1849 1891)
base 20:    66 (21 57 133 231 399 561 671 861 889 1281 1653 1729 1891 2059 2413 2501 2761 2821 2947 3059)

PhixEdit

with javascript_semantics
include mpfr.e 

function fermat_pseudoprime(integer x, base)
    if is_prime(x) then return false end if
    mpz z = mpz_init(x),
        a = mpz_init(base)
    mpz_powm_ui(z, a, x-1, z)
    return mpz_cmp_si(z,1) == 0
end function
 
constant limits = {12000, 25000, 50000, 100000}, sp = repeat('-',53)
printf(1,"Base <%s first 20 %s> <=%s\n",{sp,sp,join(limits,fmt:="%6d")})
for base=1 to 20 do
    integer count = 0, nlx = 1, nl = limits[1]
    sequence first20 = {}, counts = repeat(0,length(limits))
    for x=2 to limits[$] do
        if fermat_pseudoprime(x, base) then
            if count<20 then first20 &= x end if
            count += 1
        end if
        if x=nl then
            counts[nlx] = count
            if nlx<length(limits) then
                nlx += 1
                nl = limits[nlx]
            end if
        end if
    end for
    string s = join(first20,fmt:="%5d"),
           c = join(counts,fmt:="%6d")
    printf(1,"%2d: %s   %s\n", {base, s, c})
end for
Output:
Base <----------------------------------------------------- first 20 -----------------------------------------------------> <= 12000  25000  50000 100000
 1:     4     6     8     9    10    12    14    15    16    18    20    21    22    24    25    26    27    28    30    32    10561  22237  44866  90407
 2:   341   561   645  1105  1387  1729  1905  2047  2465  2701  2821  3277  4033  4369  4371  4681  5461  6601  7957  8321       25     38     55     78
 3:    91   121   286   671   703   949  1105  1541  1729  1891  2465  2665  2701  2821  3281  3367  3751  4961  5551  6601       25     38     53     78
 4:    15    85    91   341   435   451   561   645   703  1105  1247  1271  1387  1581  1695  1729  1891  1905  2047  2071       50     75    111    153
 5:     4   124   217   561   781  1541  1729  1891  2821  4123  5461  5611  5662  5731  6601  7449  7813  8029  8911  9881       22     35     54     73
 6:    35   185   217   301   481  1105  1111  1261  1333  1729  2465  2701  2821  3421  3565  3589  3913  4123  4495  5713       31     46     74    104
 7:     6    25   325   561   703   817  1105  1825  2101  2353  2465  3277  4525  4825  6697  8321 10225 10585 10621 11041       21     32     49     73
 8:     9    21    45    63    65   105   117   133   153   231   273   341   481   511   561   585   645   651   861   949       76    110    150    218
 9:     4     8    28    52    91   121   205   286   364   511   532   616   671   697   703   946   949  1036  1105  1288       55     83    113    164
10:     9    33    91    99   259   451   481   561   657   703   909  1233  1729  2409  2821  2981  3333  3367  4141  4187       35     53     65     90
11:    10    15    70   133   190   259   305   481   645   703   793  1105  1330  1729  2047  2257  2465  2821  4577  4921       30     41     61     89
12:    65    91   133   143   145   247   377   385   703  1045  1099  1105  1649  1729  1885  1891  2041  2233  2465  2701       35     60     91    127
13:     4     6    12    21    85   105   231   244   276   357   427   561  1099  1785  1891  2465  2806  3605  5028  5149       31     51     68     91
14:    15    39    65   195   481   561   781   793   841   985  1105  1111  1541  1891  2257  2465  2561  2665  2743  3277       33     51     69     96
15:    14   341   742   946  1477  1541  1687  1729  1891  1921  2821  3133  3277  4187  6541  6601  7471  8701  8911  9073       22     31     42     70
16:    15    51    85    91   255   341   435   451   561   595   645   703  1105  1247  1261  1271  1285  1387  1581  1687       69    102    145    200
17:     4     8     9    16    45    91   145   261   781  1111  1228  1305  1729  1885  2149  2821  3991  4005  4033  4187       31     44     63     83
18:    25    49    65    85   133   221   323   325   343   425   451   637   931  1105  1225  1369  1387  1649  1729  1921       46     69     98    134
19:     6     9    15    18    45    49   153   169   343   561   637   889   905   906  1035  1105  1629  1661  1849  1891       48     70     93    121
20:    21    57   133   231   399   561   671   861   889  1281  1653  1729  1891  2059  2413  2501  2761  2821  2947  3059       35     49     66    101

QuackeryEdit

eratosthenes and isprime are defined at Sieve of Eratosthenes#Quackery.

**mod is defined at Modular exponentiation#Quackery.

Count threshold extended to 50000.

  50000 eratosthenes

  [ over 2 < iff
      [ 2drop false ]
      done
    over isprime iff
      [ 2drop false ]
      done
    over 1 - rot **mod
    1 = ]              is fpp ( n n --> b )

  20 times
    [ say "Base: "i^ 1+ dup echo cr
      temp put
      []
      50000 times
        [ i^ 1+ temp share
          fpp if [ i^ 1+ join ] ]
      temp release
      dup size echo
      say " Fermat Pseudoprimes up to 50000"
      cr
      20 split drop
      witheach
        [ echo
          i^ 19 < if sp ]
      say "..."
      cr  cr ]
Output:
Base: 1
44866 Fermat Pseudoprimes up to 50000
4 6 8 9 10 12 14 15 16 18 20 21 22 24 25 26 27 28 30 32...

Base: 2
55 Fermat Pseudoprimes up to 50000
341 561 645 1105 1387 1729 1905 2047 2465 2701 2821 3277 4033 4369 4371 4681 5461 6601 7957 8321...

Base: 3
53 Fermat Pseudoprimes up to 50000
91 121 286 671 703 949 1105 1541 1729 1891 2465 2665 2701 2821 3281 3367 3751 4961 5551 6601...

Base: 4
111 Fermat Pseudoprimes up to 50000
15 85 91 341 435 451 561 645 703 1105 1247 1271 1387 1581 1695 1729 1891 1905 2047 2071...

Base: 5
54 Fermat Pseudoprimes up to 50000
4 124 217 561 781 1541 1729 1891 2821 4123 5461 5611 5662 5731 6601 7449 7813 8029 8911 9881...

Base: 6
74 Fermat Pseudoprimes up to 50000
35 185 217 301 481 1105 1111 1261 1333 1729 2465 2701 2821 3421 3565 3589 3913 4123 4495 5713...

Base: 7
49 Fermat Pseudoprimes up to 50000
6 25 325 561 703 817 1105 1825 2101 2353 2465 3277 4525 4825 6697 8321 10225 10585 10621 11041...

Base: 8
150 Fermat Pseudoprimes up to 50000
9 21 45 63 65 105 117 133 153 231 273 341 481 511 561 585 645 651 861 949...

Base: 9
113 Fermat Pseudoprimes up to 50000
4 8 28 52 91 121 205 286 364 511 532 616 671 697 703 946 949 1036 1105 1288...

Base: 10
65 Fermat Pseudoprimes up to 50000
9 33 91 99 259 451 481 561 657 703 909 1233 1729 2409 2821 2981 3333 3367 4141 4187...

Base: 11
61 Fermat Pseudoprimes up to 50000
10 15 70 133 190 259 305 481 645 703 793 1105 1330 1729 2047 2257 2465 2821 4577 4921...

Base: 12
91 Fermat Pseudoprimes up to 50000
65 91 133 143 145 247 377 385 703 1045 1099 1105 1649 1729 1885 1891 2041 2233 2465 2701...

Base: 13
68 Fermat Pseudoprimes up to 50000
4 6 12 21 85 105 231 244 276 357 427 561 1099 1785 1891 2465 2806 3605 5028 5149...

Base: 14
69 Fermat Pseudoprimes up to 50000
15 39 65 195 481 561 781 793 841 985 1105 1111 1541 1891 2257 2465 2561 2665 2743 3277...

Base: 15
42 Fermat Pseudoprimes up to 50000
14 341 742 946 1477 1541 1687 1729 1891 1921 2821 3133 3277 4187 6541 6601 7471 8701 8911 9073...

Base: 16
145 Fermat Pseudoprimes up to 50000
15 51 85 91 255 341 435 451 561 595 645 703 1105 1247 1261 1271 1285 1387 1581 1687...

Base: 17
63 Fermat Pseudoprimes up to 50000
4 8 9 16 45 91 145 261 781 1111 1228 1305 1729 1885 2149 2821 3991 4005 4033 4187...

Base: 18
98 Fermat Pseudoprimes up to 50000
25 49 65 85 133 221 323 325 343 425 451 637 931 1105 1225 1369 1387 1649 1729 1921...

Base: 19
93 Fermat Pseudoprimes up to 50000
6 9 15 18 45 49 153 169 343 561 637 889 905 906 1035 1105 1629 1661 1849 1891...

Base: 20
66 Fermat Pseudoprimes up to 50000
21 57 133 231 399 561 671 861 889 1281 1653 1729 1891 2059 2413 2501 2761 2821 2947 3059...

RakuEdit

use List::Divvy;
for 1..20 -> $base {
    my @pseudo = lazy (2..*).hyper.grep: { !.is-prime && (1 == expmod $base, $_ - 1, $_) }
    my $threshold = 100_000;
    say $base.fmt("Base %2d - Up to $threshold: ") ~ (+@pseudo.&upto: $threshold).fmt('%5d')
        ~ "  First 20: " ~ @pseudo[^20].gist
}
Base  1 - Up to 100000: 90407  First 20: (4 6 8 9 10 12 14 15 16 18 20 21 22 24 25 26 27 28 30 32)
Base  2 - Up to 100000:    78  First 20: (341 561 645 1105 1387 1729 1905 2047 2465 2701 2821 3277 4033 4369 4371 4681 5461 6601 7957 8321)
Base  3 - Up to 100000:    78  First 20: (91 121 286 671 703 949 1105 1541 1729 1891 2465 2665 2701 2821 3281 3367 3751 4961 5551 6601)
Base  4 - Up to 100000:   153  First 20: (15 85 91 341 435 451 561 645 703 1105 1247 1271 1387 1581 1695 1729 1891 1905 2047 2071)
Base  5 - Up to 100000:    73  First 20: (4 124 217 561 781 1541 1729 1891 2821 4123 5461 5611 5662 5731 6601 7449 7813 8029 8911 9881)
Base  6 - Up to 100000:   104  First 20: (35 185 217 301 481 1105 1111 1261 1333 1729 2465 2701 2821 3421 3565 3589 3913 4123 4495 5713)
Base  7 - Up to 100000:    73  First 20: (6 25 325 561 703 817 1105 1825 2101 2353 2465 3277 4525 4825 6697 8321 10225 10585 10621 11041)
Base  8 - Up to 100000:   218  First 20: (9 21 45 63 65 105 117 133 153 231 273 341 481 511 561 585 645 651 861 949)
Base  9 - Up to 100000:   164  First 20: (4 8 28 52 91 121 205 286 364 511 532 616 671 697 703 946 949 1036 1105 1288)
Base 10 - Up to 100000:    90  First 20: (9 33 91 99 259 451 481 561 657 703 909 1233 1729 2409 2821 2981 3333 3367 4141 4187)
Base 11 - Up to 100000:    89  First 20: (10 15 70 133 190 259 305 481 645 703 793 1105 1330 1729 2047 2257 2465 2821 4577 4921)
Base 12 - Up to 100000:   127  First 20: (65 91 133 143 145 247 377 385 703 1045 1099 1105 1649 1729 1885 1891 2041 2233 2465 2701)
Base 13 - Up to 100000:    91  First 20: (4 6 12 21 85 105 231 244 276 357 427 561 1099 1785 1891 2465 2806 3605 5028 5149)
Base 14 - Up to 100000:    96  First 20: (15 39 65 195 481 561 781 793 841 985 1105 1111 1541 1891 2257 2465 2561 2665 2743 3277)
Base 15 - Up to 100000:    70  First 20: (14 341 742 946 1477 1541 1687 1729 1891 1921 2821 3133 3277 4187 6541 6601 7471 8701 8911 9073)
Base 16 - Up to 100000:   200  First 20: (15 51 85 91 255 341 435 451 561 595 645 703 1105 1247 1261 1271 1285 1387 1581 1687)
Base 17 - Up to 100000:    83  First 20: (4 8 9 16 45 91 145 261 781 1111 1228 1305 1729 1885 2149 2821 3991 4005 4033 4187)
Base 18 - Up to 100000:   134  First 20: (25 49 65 85 133 221 323 325 343 425 451 637 931 1105 1225 1369 1387 1649 1729 1921)
Base 19 - Up to 100000:   121  First 20: (6 9 15 18 45 49 153 169 343 561 637 889 905 906 1035 1105 1629 1661 1849 1891)
Base 20 - Up to 100000:   101  First 20: (21 57 133 231 399 561 671 861 889 1281 1653 1729 1891 2059 2413 2501 2761 2821 2947 3059)

WrenEdit

Library: Wren-math
Library: Wren-gmp
Library: Wren-fmt
import "./math" for Int
import "./gmp" for Mpz
import "./fmt" for Fmt 

var one = Mpz.one

var isFermatPseudoprime = Fn.new { |x, a|
    if (Int.isPrime(x)) return false
    var bx = Mpz.from(x)
    a = Mpz.from(a)
    return a.modPow(x-1, bx) == one
}

System.print("First 20 Fermat pseudoprimes:")
for (a in 1..20) {
    var count = 0
    var x = 2
    var first20 = List.filled(20, 0)
    while (count < 20) {
        if (isFermatPseudoprime.call(x, a)) {
            first20[count] = x
            count = count + 1
        }
        x = x + 1
    }
    Fmt.print("Base $2d: $5d", a, first20)
}
var limits = [12000, 25000, 50000, 100000]
var heading = Fmt.swrite("\nCount <= $6d", limits)
System.print(heading)
System.print("-" * (heading.count - 1))
for (a in 1..20) {
    Fmt.write("Base $2d: ", a)
    var x = 2
    var count = 0
    for (limit in limits) {
        while (x <= limit) {
            if (isFermatPseudoprime.call(x, a)) count = count + 1
            x = x + 1
        }
        Fmt.write("$6d ", count)
    }
    System.print()
}
Output:
First 20 Fermat pseudoprimes:
Base  1:     4     6     8     9    10    12    14    15    16    18    20    21    22    24    25    26    27    28    30    32
Base  2:   341   561   645  1105  1387  1729  1905  2047  2465  2701  2821  3277  4033  4369  4371  4681  5461  6601  7957  8321
Base  3:    91   121   286   671   703   949  1105  1541  1729  1891  2465  2665  2701  2821  3281  3367  3751  4961  5551  6601
Base  4:    15    85    91   341   435   451   561   645   703  1105  1247  1271  1387  1581  1695  1729  1891  1905  2047  2071
Base  5:     4   124   217   561   781  1541  1729  1891  2821  4123  5461  5611  5662  5731  6601  7449  7813  8029  8911  9881
Base  6:    35   185   217   301   481  1105  1111  1261  1333  1729  2465  2701  2821  3421  3565  3589  3913  4123  4495  5713
Base  7:     6    25   325   561   703   817  1105  1825  2101  2353  2465  3277  4525  4825  6697  8321 10225 10585 10621 11041
Base  8:     9    21    45    63    65   105   117   133   153   231   273   341   481   511   561   585   645   651   861   949
Base  9:     4     8    28    52    91   121   205   286   364   511   532   616   671   697   703   946   949  1036  1105  1288
Base 10:     9    33    91    99   259   451   481   561   657   703   909  1233  1729  2409  2821  2981  3333  3367  4141  4187
Base 11:    10    15    70   133   190   259   305   481   645   703   793  1105  1330  1729  2047  2257  2465  2821  4577  4921
Base 12:    65    91   133   143   145   247   377   385   703  1045  1099  1105  1649  1729  1885  1891  2041  2233  2465  2701
Base 13:     4     6    12    21    85   105   231   244   276   357   427   561  1099  1785  1891  2465  2806  3605  5028  5149
Base 14:    15    39    65   195   481   561   781   793   841   985  1105  1111  1541  1891  2257  2465  2561  2665  2743  3277
Base 15:    14   341   742   946  1477  1541  1687  1729  1891  1921  2821  3133  3277  4187  6541  6601  7471  8701  8911  9073
Base 16:    15    51    85    91   255   341   435   451   561   595   645   703  1105  1247  1261  1271  1285  1387  1581  1687
Base 17:     4     8     9    16    45    91   145   261   781  1111  1228  1305  1729  1885  2149  2821  3991  4005  4033  4187
Base 18:    25    49    65    85   133   221   323   325   343   425   451   637   931  1105  1225  1369  1387  1649  1729  1921
Base 19:     6     9    15    18    45    49   153   169   343   561   637   889   905   906  1035  1105  1629  1661  1849  1891
Base 20:    21    57   133   231   399   561   671   861   889  1281  1653  1729  1891  2059  2413  2501  2761  2821  2947  3059

Count <=  12000  25000  50000 100000
------------------------------------
Base  1:  10561  22237  44866  90407 
Base  2:     25     38     55     78 
Base  3:     25     38     53     78 
Base  4:     50     75    111    153 
Base  5:     22     35     54     73 
Base  6:     31     46     74    104 
Base  7:     21     32     49     73 
Base  8:     76    110    150    218 
Base  9:     55     83    113    164 
Base 10:     35     53     65     90 
Base 11:     30     41     61     89 
Base 12:     35     60     91    127 
Base 13:     31     51     68     91 
Base 14:     33     51     69     96 
Base 15:     22     31     42     70 
Base 16:     69    102    145    200 
Base 17:     31     44     63     83 
Base 18:     46     69     98    134 
Base 19:     48     70     93    121 
Base 20:     35     49     66    101