6

A Carmichael number is a composite number $N$, such that $a^{N-1}\equiv 1\mod N$ holds for every $a$ coprime to $N$. $N$ is a Carmichael number if

  • $N$ is odd and squarefree
  • $N$ has at least three distinct prime factors
  • For each prime factor $p|N$ we have $p-1|N-1$

A positive integer $M$ is given.

How can I efficiently find the smallest Carmichael number $N\ge M$ ?

For example, the smallest Carmichael-number above $10^{10}$ is $$10017089857=73\cdot163\cdot577\cdot1459$$

I am looking for a method more efficient than brute force, something like a Carmichael-sieve. Any ideas ?

Peter
  • 86,576
  • 1
  • The Feitsma/Galway list (http://www.cecm.sfu.ca/Pseudoprimes/index-2-to-64.html) has all results to 2^64. You can pull them out of the annotated file, or do a is_carmichael test on the psps. BTW, a brute force method for me takes about 2.5s to find next-carmichael(10^10), but has bad growth (e.g. 4s for 10^11, 25s for 10^12, 140s for 10^13). – DanaJ May 26 '17 at 18:31
  • @DanaJ The best approach I found yet is to start checking whether a number is a Fermat-pseudoprime to base $2$. Is there any better method ? – Peter May 30 '17 at 08:04
  • 1
    @Peter for me: 1: reject if < 561 or even, 2: reject if divisible by 9,25,49,121,169, 3: reject using Korselt's criterion for small divisors (e.g. reject if n div 5 and n-1 not div 4, reject if n div 7 and n-1 not div 6, etc.), 4: reject if base 2 pseudoprime, 5: the long way... factor, reject if num factors < 3, reject if any factor duplicated, check Korselt for each factor. – DanaJ May 30 '17 at 16:06

2 Answers2

9

I don't believe that there is a much better way to efficiently find a specific Carmichael number. Actually, the computations done at OEIS by Pinch list all Carmichael numbers up to $1713045574801$. So $10^{10}$ is not a problem.
On the other hand one can find an easy upper bound for the smallest Carmichael number above $10^{10}$ by using Chernick's criterion: if $k$ is a positive integer such that $6k + 1, 12k + 1$, and $18k + 1$ are all prime then the product $$ n = (6k + 1)(12k + 1)(18k + 1) $$ is a Carmichael number. For $k=195$ the criterion applies, because $6\cdot 195+1$, $12\cdot 195+1$, and $18\cdot 195+1$ are all prime, so that $f(195)=9624742921$ is a lower bound for the largest Carmichael number below $10^{10}$ (which is $9999109081$). For a Carmichael number above $10^{10}$ we can take $k=206$, so that $11346205609$ is a Carmichael number above $10^{10}$. Of course, this will not be enough to efficiently find a specific Carmichael number as said.

Dietrich Burde
  • 140,055
  • 1
    This is not my meat, but it’s a wonderful answer. – Lubin May 24 '17 at 19:11
  • A Sophie Germain prime is a prime p such that 2p-1 is prime. If 6k+1 and 12k+1 are primes then 6k+1 is a Sophie Germain prime. It is not known whether there is a largest one. If there is, then this method has finite scope. I do not know the proof that there is no largest Carmichael number. Perhaps the proof has something in it on this Q. – DanielWainfleet May 24 '17 at 20:06
  • 1
    @DanielWainfleet: Sophie Germain primes are $2p\color{red}{+}1$. – ccorn May 25 '17 at 02:15
  • @Lubin The method gives a reasonable upper bound (assuming that there are infinite many $k$, such that $6k+1,12k+1,18k+1$ are prime, which is not proven , as far as I know, but very likely true). But for large $M$, the bound is not tight enough to finish the calculation with brute force. – Peter May 25 '17 at 09:08
  • @DanielWainfleet I do not think that the proof that infinite many Carmichael numbers exist can be derived from this criterion. It is likely, but as far as I know not proven that the above formula produces infinite many Carmichael numbers. – Peter May 25 '17 at 09:10
1

Here's an algorithm in Python code for efficiently computing the Carmichael numbers in a given range:

from sympy.ntheory import sieve, isprime, prime
from sympy.core import integer_nthroot
from math import lcm, gcd, isqrt

def find_carmichael_numbers_in_range(A, B): max_p = 1+((1 + isqrt(8*B + 1)) >> 2)

def _rec(m, lam, lo, k):

    if k == 1:

        lo = max(lo, A // m + (1 if A % m else 0))
        hi = min(B // m + 1, max_p)

        u = pow(m, -1, lam)
        j,r = divmod(lo - u, lam)
        if r: j += 1
        u += j*lam

        for p in range(u, hi, lam):
            if (m*p - 1) % (p - 1) == 0 and isprime(p):
                yield m*p

    else:
        hi = (integer_nthroot(B // m, k))[0]+1
        for p in sieve.primerange(lo, hi):
            if gcd(m, p - 1) == 1:
                yield from _rec(m*p, lcm(lam, p - 1), p + 2, k - 1)

def _rec1():
    k = 3
    l = 3*5*7
    while l &lt; B:
        yield from _rec(1, 1, 3, k)
        k += 1
        l *= prime(k+1)

return sorted(_rec1())

assert find_carmichael_numbers_in_range(561, 2000) == [561, 1105, 1729] assert find_carmichael_numbers_in_range(562, 2000) == [1105, 1729] assert find_carmichael_numbers_in_range(10585, 75361) == [10585, 15841, 29341, 41041, 46657, 52633, 62745, 63973, 75361] assert [len(find_carmichael_numbers_in_range(1, 10**n)) for n in range(1, 8)] == [0, 0, 1, 7, 16, 43, 105] # A055553

print(find_carmichael_numbers_in_range(1010, 1010 + 10**8)) #=> [10017089857, 10017782401, 10028287081, 10031651841, 10054063041, 10073377801, 10087556401]

Trizen
  • 101