# Novel Prime Factors of Repunits

### The Background

Folks were talking prime factorization in chat and we found ourselves talking about repunits. Repunits are a subset of the numbers known as repdigits, which are numbers consisting of only repeating digits, like 222 or 4444444444444444, but repunits consist only of 1.

The first couple repunits are therefore 1, 11, 111, etc. These are referred to by Rn, so R1=1, R2=11, etc., and are generated by the formula R(n) = (10^n - 1)/9, with n > 0.

Prime factorization of these repunit numbers follows sequence A102380 in the OEIS. For example:

R1 = 1
R2 = 11
R3 = 111 = 3 * 37
R4 = 1111 = 11 * 101
R5 = 11111 = 41 * 271
R6 = 111111 = 3 * 7 * 11 * 13 * 37
R7 = 1111111 = 239 * 4649
...

### The Challenge

Write a program or function which, when given an input integer n with n >= 2 via STDIN or equivalent, outputs or returns the novel prime factors for Rn, in any convenient format. "Novel prime factors" here means all x where x is a prime factor of Rn, but x is not a prime factor for any previous Rk, with 1 <= k < n (i.e., if we write the prime factors for all R in sequence, we've not seen x before).

### The Examples

Input: 6
Output: 7, 13
(because 3, 11, and 37 are factors of a smaller R_k)

Input: 19
Output: 1111111111111111111
(because R_19 is prime, so no other factors)

Input: 24
Output: 99990001
(because 3, 7, 11, 13, 37, 73, 101, 137, 9901 are factors of a smaller R_k)

Input: 29
Output: 3191, 16763, 43037, 62003, 77843839397
(because no factors of R_29 are also factors of a smaller R_k)


### The Extras:

• Your code can do anything or nothing if n < 2.
• You can assume a "reasonable" upper limit for n for testing and execution purposes -- your code will not be expected to output for n = 10000000, for example, but your algorithm should work for such a case if given unlimited computing power and time.
• Here is a site dedicated to the factorizations of repunits for reference.
• I've not gone through the math, but I propose a hypothesis that every n has a distinct result for this algorithm -- that is, no n exists such that Rn has no novel factors. I'll offer a 250-point bounty if someone proves or disproves that in their answer. Thomas Kwa proposed an elegant proof, and I awarded the bounty.
• Built-in prime checking and factorisation are fair game? – Martin Ender Dec 14 '15 at 16:21
• @MartinBüttner No restrictions. – AdmBorkBork Dec 14 '15 at 16:21
• OEIS A204845 – alephalpha Dec 15 '15 at 4:39
• @alephalpha Because of course there's an OEIS link ;-) Thanks! – AdmBorkBork Dec 15 '15 at 13:55

# Pyth, 1614 13 bytes

-F_m{P/^Td9SQ


As best I can tell, 19 takes forever.

-F_m{P/^Td9SQ      Implicit: Q = input
SQ      Inclusive range 1 to Q
Implicit lambda d:
{P                  unique primes dividing
^Td              10**d
/   9                  //9
m               map lambda over the range
m{P/^Td9SQ      Unique prime factors of all repunits up to the Qth
_                Reverse the list
-F                 Reduce by set difference


Try it here.

• If you run it from the command line and have SymPy installed, 19 completes in under a second. The first one that takes over 2 seconds is the input 38. – isaacg Dec 14 '15 at 19:17

# Proof that every repunit has a novel prime factor

Using Zsigmondy's Theorem, the proof is simple. From Wikipedia:

In number theory, Zsigmondy's theorem, named after Karl Zsigmondy, states that if a > b > 0 are coprime integers, then for any integer n ≥ 1, there is a prime number p (called a primitive prime divisor) that divides an − bn and does not divide ak − bk for any positive integer k < n, with the following exceptions: [things that don't apply here].

Consider the repunits times 9: that is, 10n-1. By Zsigmondy's Theorem with a=10, b=1, there is some prime p | 10n-1 that does not divide any 10k-1, k < n.

• Since p is prime and 10n-1 = 9 · Rn, it must divide either 9 or Rn.

• p cannot divide 9, since 9 = 101-1.

• Therefore p divides Rn.

• p cannot divide any Rk, since it doesn't divide 10k-1 = 9 · Rk.

Therefore the p from Zsigmondy's Theorem is a novel prime factor of any Rn, n ≥ 2. ∎

## An example of a repeated novel prime factor

The prime 487 is a repeated prime factor of R486:

By the Fermat-Euler theorem, 10487-1 ≡ 1 (mod 487). The order of 10 mod 487, that is, the smallest power of 10 that is 1 mod 487, therefore must be a divisor of 486. In fact, the order is equal to 486. It also happens that not only is 10486 ≡ 1 (mod 487), it's also 1 (mod 4872).

See here that the order of 10 mod 487 is 486; that is, that no smaller 10k-1 is divisible by 487. Obviously, 487 doesn't divide 9, so it must divide R486.

## CJam, 18 bytes

{{)'1*imf}%W%:-_&}


Test it here.

Starting at 19 (the first repunit prime after 2) it will take a fairly long time.

### Explanation

This is an unnamed block (CJam's equivalent of a function), which expects the input number n on the stack and leaves a list of prime factors instead:

{      e# Map this block over [0 1 ... n-1]...
)'1* e#   Increment and create a string of that many 1s.
i    e#   Convert to integer.
mf   e#   Get its prime factors.
}%
W%     e# Reverse the list.
:-     e# Fold set difference onto the list, removing from the first list the elements of
e# all other lists.
_&     e# Remove duplicates. Unfortunately, this necessary. Thomas Kwa found that the 486th
e# repunit contains 487^2 (where 487 is a novel prime factor).

• Did ... did you seriously just golf from 20 to 16 in the past three minutes? >.> – AdmBorkBork Dec 14 '15 at 16:31
• @TimmyD Sort of... I had to go up to 18 again for now, because it turned out my code was based on an assumption which I can neither prove or disprove right now. – Martin Ender Dec 14 '15 at 16:42
• Ooo, that's an interesting case - duplicate novel factors. Nice catch. – AdmBorkBork Dec 14 '15 at 17:05

import Data.Numbers.Primes
f n=[x|x<-primeFactors$div(10^n-1)9,notElem x$f=<<[1..n-1]]


Usage example: f 8 -> [73,137].

Takes a lot of time and memory for big n.

Direct implementation of the definition: take all prime factors x of Rn which do not show up before (f=<<[1..n-1] are all prime factors of R1 ... R(n-1)).

# Mathematica 82 74 63 bytes

With 11 bytes saved thanks to alephalpha.

Complement@@Reverse@FactorInteger[(10^Range@#-1)/9][[;;,;;,1]]&


Prime factors of R70

(10^70 - 1)/9 = 1111111111111111111111111111111111111111111111111111111111111111111111

FactorInteger[(10^70 - 1)/9]


{{11, 1}, {41, 1}, {71, 1}, {239, 1}, {271, 1}, {4649, 1}, {9091, 1}, {123551, 1}, {909091, 1}, {4147571, 1}, {102598800232111471, 1}, {265212793249617641, 1}}

Novel prime factors of R70

Complement@@Reverse@FactorInteger[(10^Range@#-1)/9][[;;,;;,1]]&[70]


{4147571, 265212793249617641}

• Complement@@Reverse@FactorInteger[(10^Range@#-1)/9][[;;,;;,1]]& – alephalpha Dec 15 '15 at 4:20
• Kindly explain the meaning of [[;;,;;,1]] or [[1 ;; All, 1 ;; All, 1]]. I'm puzzled! – DavidC Dec 15 '15 at 13:22
• @DavidCarraher It takes the first element of every element of every element of a list. – LegionMammal978 Dec 15 '15 at 13:25
• @DavidCarraher [[;;,;;,1]] is the same as [[All,All,1]]. – alephalpha Dec 15 '15 at 13:46
• Now that makes sense. BTW, your re-location of Range was very clever. – DavidC Dec 15 '15 at 20:29

# MATL, 25 bytes

This works for input up to 16:

10,i:^9/Y[t0)Yftb!w\~s1=)


The following version uses 31 bytes and works up to 18. For 19 it requires about 4 GB memory (I haven't been able to run it).

10,i:^9/Y[t0)5X2Y%Yfotb!w\~s1=)


### Example

>> matl
> 10,i:^1-,9/t0)5X2Y%Yfotb!w\~s1=)
>
> 6
7 13


### Explanation

Consider for concreteness input 6. First the prime divisors of 111111 are computed; in this case the results are 3, 7, 11, 13, 37. Then the modulo operation (division remainder) is computed for all combinations of numbers 1, 11, ... 111111 and the computed divisors. This exploits MATL's implicit singleton expansion. The result is in this case a 6x5 matrix, with each column corresponding to one of the divisors. Accepted divisors (columns) are those for which only 1 value (namely the last) gives zero remainder.

10,i:^9/Y[   % generate vector with 1, 11, ... depending on input number, say "n"
t0)          % pick the last element: 111...1 (n ones)
5X2Y%        % * convert to uint64, so that larger numbers can be handled
Yf           % prime factors
o            % * convert to double precision, so that modulus can be done
t            % duplicate
b            % bubble up element in stack
!            % transpose
w            % swap elements in stack
\            % modulus after division (element-wise, singleton expansion)
~s           % number of zero values in each column
1=           % is equal to 1? (element-wise, singleton expansion)
)            % index divisors with that logical index


(*) Removed in short version

• That's a clever way of doing it. – AdmBorkBork Dec 14 '15 at 21:55

# Julia, 103 bytes

R(n)=[keys(factor((10^n-1)÷9))...]
n->(r=R(n);isprime(r)?r:setdiff(r,reduce(vcat,[R(i)for i=1:n-1])))


This is an unnamed function that calls a helper function R. To call it, give the main function a name, e.g. f=n->....

Ungolfed:

function R(n::Integer)
collect(keys(factor((10^n - 1) ÷ 9)))
end

function f(n::Integer)
r = R(n)
if isprime(r)
r
else
setdiff(r, reduce(vcat, [R(i) for i = 1:n-1]))
end
end


# LabVIEW, 33 LabVIEW Primitives

19 takes forever...

Work by saving all Primes and deleting elements from the last set when they are found in the other array.

# J, 24 bytes

[:({:-.}:)@:q:[:+/\10^i.


Expects extended-precision numerics after 6 (e.g. 19x instead of 19).

Try it online!

There's probably a shorter way to generate the repunits which avoids the caps, too.

# Explanation

[: ({: -. }:) @: q: [: +/\ 10 ^ i.
i. Range [0, input)
10 ^    10 raised to the power of the range
+/\         Running sum of this list (list of repunits)
q:                Prime factors of the repunits
@:                   Composed with
({: -. }:)                      Unique prime factors of last repunit
{:                              Factors of last repunit (tail of matrix)
-.                           Set subtraction with
}:                        Rest of the matrix (curtail)


### How it works, visually

I think that these sorts of visual explanations are easier to stomach for those who don't know J. These are results from the REPL.

   10 ^ i.6
1 10 100 1000 10000 100000

+/\ 10 ^ i. 6
1 11 111 1111 11111 111111

q: +/\ 10 ^ i. 6
0   0  0  0  0
11   0  0  0  0
3  37  0  0  0
11 101  0  0  0
41 271  0  0  0
3   7 11 13 37

({: -. }:) q: +/\ 10 ^ i. 6
7 13