Factorization of Sylvester's sequence

Sylvester's sequence is an integer sequence (si) where s0 is 2, and si is one plus the product of all i previous terms: si = s0s1...si-1 + 1.
The sequence begins 2, 3, 7, 43, 1807, 3263443, 10650056950807, 113423713055421844361000443.
The empty product of no terms is 1, so s0 = 1 + 1 = 2 is a natural start.
It can easily be proved by induction that si+1 = si (si−1) + 1. This reminds of the Fermat numbers Fi = 22i+1, where:
     F0 = 3,    Fi = F0F1...Fi-1 + 2,   Fi+1 = (Fi−1) (Fi−1) + 1.
The numbers in Sylvester's sequence are relatively prime. Most primes do not divide any of them. It is not known whether all the numbers are square-free. No prime factor at this site has a square that also divides the number.

The following table shows known factors up to s22. Pn is an n-digit prime. Cn is an n-digit composite number with unknown factorization.

i

Factorization of si

0 2 (prime)
1 3 (prime)
2 7 (prime)
3 43 (prime)
4 13 139
5 3263443 (prime)
6 547 607 1033 31051
7 29881 67003 9119521 6212157481
8 5295435634831 31401519357481261 77366930214021991992277
9 181 1987 112374829138729 114152531605972711 P68
10 2287 2271427 21430986826194127130578627950810640891005487 P156
11 73 C416
12 2589377038614498251653 2872413602289671035947763837 C785
13 52387 5020387 5783021473 401472621488821859737 287001545675964617409598279 C1600
14 13999 74203 9638659 57218683 10861631274478494529 C3293
15 17881 97822786011310111 54062008753544850522999875710411 C6618
16 128551 C13335
17 635263 1286773 21269959 C26661
18 50201023123 139263586549 60466397701555612333765567 C53313
19 775608719589345260583891023073879169 C106685
20 352867 6210298470888313 C213419
21 387347773 1620516511 C426863
22 91798039513 C853750

s0 to s9 are easy to factor. Ken Takusagawa made the largest complete factorization by finding a 44-digit factor of s10 with GMP-ECM.
Sean A. Irvine found several large factors from s12 to s18 and added them to Wikipedia. He has reported the following GMP-ECM work to me. s11: 4590 curves at B1=11e6. s12: 1000 curves at 1e6. s13: 100 at 50k. s14: 100 at 10k. s15: 100 at 10k. s16: 100 at 10k. Only few curves above s16.
Dario Alpern found the two known factors of s21 in 2006.
I trial factored to 1012 in 2007 which gave the new factor 50201023123 in s18 (Irvine had already found a larger factor), and the first known factor of s22.
Christophe Clavier found the factor 6210298470888313 of s20 in 2007, the factor 54062008753544850522999875710411 of s15 in 2009, the factor 60466397701555612333765567 of s18 in 2012, and the factor 775608719589345260583891023073879169 of s19 in 2012.
The smallest number without a known factor is s30 which has 218562698 digits and unknown primality status.
The large composite divisors were proven composite with PrimeForm/GW. A base 3 Fermat prp test was performed with these results:

s18: 1080204337.....2994420807/(50201023123*139263586549*60466397701555612333765567) is composite: RES64: [439B8ADB88D2121B] (200.7818s+0.0867s)
s19: 1166841411.....6400110443/775608719589345260583891023073879169 is composite: RES64: [A3CF26A0FD91EDAC] (848.6814s+25.3268s)
s20: 1361518878.....6197545807/(352867*6210298470888313) is composite: [2FAA731D98F69934] (3957.6467s+0.4931s)
s21: 1853733657.....3665735443/(387347773*1620516511) is composite: [1ED1B20FEACEAE0B] (17885.7170s+1.9189s)
s22: 3436328472.....4400670807/91798039513 is composite: [1A69A90C683A93A6] (83768.4478s+8.6656s)
The below list shows known factors of s23 to s200. The unfactored part has unknown primality status in all cases. According to heuristics based on the fast growth, it is unlikely that any si above s5 is prime.
In 1991 Ilan Vardi listed the prime factors below 5107 up to s200. He also showed (I don't know whether he was first) that all prime factors above 3 are of form 6n+1. Recall si+1 = si (si−1) + 1 and consider what happens modulo p when going from si to si+1. If n(n−1)+1 ≡ 0 (mod p) has no solution, then p cannot be a factor of si+1. This happens when there is no solution to −3 ≡ n2 (mod p), and it has no solution for a prime p iff p is of form 6n−1.
Next, consider si+2 = si+1 (si+1−1) + 1 = (si(si−1)+1) (si(si−1)) + 1. If (n(n−1) +1) (n(n−1)) + 1 ≡ 0 (mod p) has no solution, then p cannot be a factor of si+2. This reasoning can be continued with si+3 and so on, but the ratio of potential prime factors which is eliminated will decrease each time.
Dario Alpern computed the prime factors below 2109 up to s100 in 2006. I have not seen Vardi's list and had not seen Alpern's when I made my own trial factoring in 2007 to 1012 for s23 to s25, and to 1011 for s26 to s200. Christophe Clavier found the factor 206866598749591633 of s24 in 2012.
 
 i  Known prime factors of s_i
--- --------------------------
 23 74587
 24 206866598749591633
 25 624116543851
 26 27061
 27 164299, 3229081, 16194353551
 28 20929, 2621897371
 29 1171, 298483, 97562299
 30
 31 1679143
 32
 33
 34 120823, 447841
 35 2408563
 36 38218903, 818864029
 37 333457, 117602053
 38 30241
 39 4219, 1372512301
 40 1085443, 156507607, 171818923
 41 7603
 42 1861, 84819199
 43
 44 23773, 290791
 45 51769
 46 1285540933
 47 429547, 196973983
 48
 49 8323570543
 50
 51 179585227
 52
 53 9829, 15238744627
 54
 55
 56 47101199947
 57 763176709, 15517465963
 58
 59 2029, 38329, 320329, 3567469
 60
 61 50023, 1445860807
 62
 63
 64 1459, 11234941
 65 72091, 609421, 6377833117
 66
 67 245563, 3346633, 343870981
 68 6367063, 10516273813
 69 12763
 70 384061
 71 3799489, 8401963
 72
 73 11844787
 74 35869, 26664046099
 75 20161, 42428887
 76 123427, 160200259
 77
 78
 79 51973, 195502537
 80 2437, 49848925897
 81
 82 6230971, 57741984199
 83
 84 24049307239
 85
 86 9436201
 87
 88 151477
 89 41077, 2171593
 90
 91
 92 15373, 21661, 1386013
 93 429446209
 94 12593197
 95 535879
 96 1407223, 1143074323
 97 2748116731
 98 20479
 99
100 537037, 142997167
101 73381873
102
103 4519
104 13291
105 9123469
106 97549
107 1101733, 116289529
108 297589, 549733
109 118953139
110 943567
111 1285633
112
113 132637, 45352921
114 473611
115 6906223
116 148207, 51205591, 331215169
117 17977
118
119 54919, 644702497
120
121
122
123 1944721, 234341827
124
125
126 32779, 1604059579
127
128
129
130 26778949, 142670533
131 36680803
132
133
134 2541461581
135 256204723, 315917071
136 27031
137 8221, 149197, 186749617
138 769858669, 11304677191
139 46380739
140 6469
141 6886081
142 1182253, 2211019
143
144 615721
145 2564857, 3113936041
146 13147, 8186719
147
148 7205893
149
150
151 2686393, 369582877
152
153
154
155
156 232051, 438049093
157
158
159
160
161 64621, 1530229, 8244391
162
163
164
165
166 1763911, 6051648259
167
168
169 695131, 3503498041
170 215899
171
172 35099257
173 839437, 792599053
174
175 46022962369
176 31159, 28384141
177
178 4508756437
179
180 57853
181
182
183 107053
184 60919, 505123
185 39360469
186
187
188 7948981
189
190 2986572571
191
192
193 12436189
194 80173, 11073379, 734468401
195 1156390513, 61625866987
196
197 562909, 1914031477
198 13163791, 24518173
199 19597, 150649, 95037412693
200 765241

sylvester-factors.txt contains all known prime factors, written as [i, p] where p divides si, and sorted by p.
The factors include: All factors below 232 for all i (this page stops at i = 200), all factors below 1011 for i ≤ 200, all factors below 1012 for i ≤ 25, all factors of any size for i ≤ 10, and some factors above 1012 for 11 ≤ i ≤ 24.

Let p be a fixed prime. The recursive formula si+1 = si (si−1) + 1 means that the sequence (si) modulo p must eventually enter a cycle since there is a finite number of possible values modulo p. Therefore, if a cycle in (si) modulo p is detected before reaching a number divisible by p, then p will never be a factor in Sylvester's sequence.
Vardi reported that 1166 out of the first three million primes (meaning up to 49979687 with no further factors below 5107) are prime factors of some si. My triple-checked computation says 1167 factors (the first 1167 entries in sylvester-factors.txt), corresponding to 1 in 2571. My additional computations say that 8181 of the 203280221 primes below 232 are factors (1 in 24848).

Here is a PARI/GP function to test whether p divides si, using modular multiplication without computing si:

isfactor(i,p) = local(j,s);s=2%p;for(j=1,i,s=(s*(s-1)+1)%p);s==0
It is slow and unsuited to search new factors but it can verify all known factors in 20 GHz minutes with this PARI/GP function:
verify_factors(file) = {
  local(f,c,j,p,n);
  f=readvec(file);
  c=j=p=0;
  for(n=1,#f,
    if(isfactor(f[n][1],f[n][2]), c++);
    if(isfactor(f[n][1],(f[n][2])^2), j++);
    if(isprime(f[n][2]), p++);
    if (n>1 && f[n][2]<=f[n-1][2], print("Number "n" not increasing"));
  );
  print(#f" tested, "p" primes, "c" correct factors, "#f-c" incorrect, "j" square factors.");
}

? verify_factors("sylvester-factors.txt")
8223 tested, 8223 primes, 8223 correct factors, 0 incorrect, 0 square factors.

Links
Wikipedia, Sylvester's sequence.
MathWorld by Eric W. Weisstein, Sylvester's Sequence.
Ken's blog by Ken Takusagawa, Factoring Sylvester's sequence and Sylvester 10th factored.
Mersenneforum.org, Number 100 ("alpertron" is Dario Alpern) and Sylvester's Sequence.

Reference copied from MathWorld (without having access to it):
Vardi, I. "Are All Euclid Numbers Squarefree?" and "PowerMod to the Rescue." §5.1 and 5.2 in Computational Recreations in Mathematica. Reading, MA: Addison-Wesley, pp. 82-89, 1991.

This page is at http://users.cybercity.dk/~dsl522332/math/sylvester-factors.htm and licensed under the GFDL.
Page by Jens Kruse Andersen, jens.k.a@get2net.dk   home
First uploaded 27 August 2007. Last updated 9 September 2012