## Lucas-Carmichael Numbers

### January 6, 2015

We start the new year with a simple task from number theory. A Lucas-Carmichael number is a positive composite integer *n*, odd and square-free, such that *p* + 1 divides *n* + 1 for all prime factors *p* of *n*. For instance, 2015 is a Lucas-Carmichael number because 2015 = 5 × 13 × 31 and 5 + 1 = 6, 13 + 1 = 14, and 31 + 2 = 32 all divide 2015 + 1 = 2016. The restriction that *n* be square-free prevents the cubes of prime numbers, such as 8 or 27, from being considered as Lucas-Carmichael numbers.

Your task is to write a program to find all the Lucas-Carmichael numbers less than a million. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Perl, 5s for 1e7 (135 L-C numbers), <1.5min for 1e8 (323 L-C numbers). Perl:

use warnings;

use strict;

use ntheory qw/:all/;

foroddcomposites {

my $n=$_;

say unless moebius($n) == 0 || vecsum( map { ($n+1) % ($_+1) != 0 } factor($n) );

} 1e6;

I decided to use a sieve.

For p = 3, write out the multiples of p and of (p+1):

3 6 9 12 15 18 21 24 27

4 8 12 16 20 24 28

The only values of n for which (n+1)%(p+1) == 0 are n = 15, 27, 39, …. Similarly, for p = 5, the only values are n = 35, 65, 95, …. In general, for a prime ‘p’, the only values of ‘n’ for which (n+1)%(p+1) == 0 are given by: p + k*p*(p+1) for k = 1,2,3,…

Initialize the sieve to all 1s.

For each odd prime ‘p’:

for k in 1, 2, 3, …:

sieve[p+k*p*(p+1)] *= p

return a list of all n where sieve[n] == n

The code below sieves the primes I parallel with sieving the LCNs.

The sieves also only contain odd n, so the index for a value of n is given by n//2.

Finds 323 LCNs under 1e8 in about 12 seconds.

Haskell

import Data.List

factors :: Int -> [Int]

factors n = factors’ n 2

where

factors’ n f

| f*f > n = [n]

| n `mod` f == 0 = f : factors’ (n `div` f) f

| otherwise = factors’ n (f+1)

noDupes :: [Int] -> Bool

noDupes xs = length(nub xs) == length xs

lcn :: [Int] -> [Int]

lcn xs = filter (\x ->

let fs = factors x

in

length(fs) > 1 && noDupes(fs) && (all (\f -> (x + 1) `mod` (f + 1) == 0) fs) ) xs

limit = 1000000

n = [3,5..limit]

main = print (lcn n)