Multiplicative prime number sets

[henryzz]

I had the idea for generating sets of primes this morning:

1. Start with the primes 2 and p.
2. Multiply the numbers in the set(2 and p in first iteration) and find all divisors(including composite divisors).
3. Subtract 1 from all these proper divisors
4. Add any new primes to the set of primes and go back to step 2.

Starting with 2 and 7, I ended up with the set {2 7 13 181 32941 11924641 edit: missed prime 65881. There are then many more terms which I haven't fully calculated yet.}. 2 and 19 results in {2 19 37 73}

Is there a limit to the size of these sets? The number of new potential primes increases every time a new prime is added, but the probability of primality also decreases as the numbers generally increase in size.

It is possible to replace step with addition by 1 rather than subtraction.

[Dr Sardonicus]

If widely-believed conjectures are true, there is no upper bound on the number of primes that can be produced.

Note that if p is an odd prime, the divisors of 2*p are 1, 2, p, and 2*p. Subtracting 1 gives 0, 1, p-1 and 2*p - 1. So you are dead in the water unless p and 2*p - 1 are both prime. "Everybody knows" that there are infinitely many such p, but AFAIK nobody knows how to prove it. There is also a conjectured asymptotic formula.

The next iteration then produces 4*p - 3. Again, it is conjectured that there are infinitely many p such that p, 2*p - 1, and 4*p - 3 are all prime, with a conjectured asymptotic formula.

This can be continued indefinitely. For any positive integer n, there are (conjecturally) infinitely many p for which p and 2^k*p + 2^k - 1 are all prime, for k = 1 to n.

The asymptotic formula picks up an extra log factor in the denominator for each added term.

Of course, for any given p, the sequence of consecutive primes will be finite (although I'm too lazy at the moment to try to prove this).

The longest such sequences I found in a desultory search were with p = 15514861 and p = 62573941; these continue until 128*p - 127.

EDIT: The sort of sequence I gave is well known. It is called a Cunningham chain of the second kind. (The alternate problem where you add 1 leads to "Cunningham chains of the first kind.")

[henryzz]

I should have thought of Cunningham Chains. That leaves the remaining question do most sets grow to very large sizes or run out of new primes?

I need to upgrade my code for this to use previous primes rather than factor the large composite each time which quickly becomes unsustainable.

fordiv(2*7*13*181*32941*65881*922333*11924641*2337191821*394973855821*850699085221*5922444074413,x,if(isprime(x-1),print(x-1)))

[Dr Sardonicus]

Of course, most sets run out of primes immediately. That is, for most primes p, 2p - 1 is composite. Starting with 2 and 7, I ran into an embarrassment of riches:

? s=Set([2,7]);for(k=1,6,m=#s;v=vecsort(eval(s));print(k" "#v" "v);print();u=vector(m,i,1)~;M=Mat([v~,u]);fordiv(M,d,p=d-1;if(ispseudoprime(p),s=setunion(s,Set(p)))));v=vecsort(eval(s));print(7" "#v" "v)

1 2 [2, 7]

2 3 [2, 7, 13]

3 4 [2, 7, 13, 181]

4 5 [2, 7, 13, 181, 32941]

5 7 [2, 7, 13, 181, 32941, 65881, 11924641]

6 10 [2, 7, 13, 181, 32941, 65881, 922333, 11924641, 394973855821, 284389833087001]

7 29 [2, 7, 13, 181, 32941, 65881, 922333, 11924641, 2337191821, 394973855821, 850699085221, 71885241759421, 142980380787217, 142980535807201, 284389833087001, 9419842870102370521, 1714383664924632549457, 3391246948002241778641, 1572571685029287436525959493, 20331412254982256557465004941, 43789360924492372446131259181, 1860295829820854184468425439541, 7359859521855984153835211751181, 121890303361521579271853639370932461, 487553284314997301558832816949395457, 88732004080775282379279070846022297137, 8647957628487867905734352520146942344141, 31944528501630323994538127586861144992641, 2906952093648359483502969610404364194330421]
?

Uh-oh!

EDIT: I devised a way to push a bit further. In theory, one could set the upper bound on the loop counter i to be 2^29 - 1 rather than 2^28 + 100000, and get the complete set of primes from the eighth iteration. In practice, it might take quite a while. The ninth iteration would appear to be computationally beyond reach.

? v=[2, 7, 13, 181, 32941, 65881, 922333, 11924641, 2337191821, 394973855821, 850699085221, 71885241759421, 142980380787217, 142980535807201, 284389833087001, 9419842870102370521, 1714383664924632549457, 3391246948002241778641, 1572571685029287436525959493, 20331412254982256557465004941, 43789360924492372446131259181, 1860295829820854184468425439541, 7359859521855984153835211751181, 121890303361521579271853639370932461, 487553284314997301558832816949395457, 88732004080775282379279070846022297137, 8647957628487867905734352520146942344141, 31944528501630323994538127586861144992641, 2906952093648359483502969610404364194330421];s=Set(v);for(i=2^28+1,2^28+100000,w=binary(i);f=1;for(j=1,29,if(w[j]==1,f*=v[j]));p=f-1;if(ispseudoprime(p),s=setunion(s,Set(p))))
? print(#s)
610
? v=vecsort(eval(s));
? print(v[610])
1271263558764320475807623718922557314540894891927785434354432334853043584019534069482685269646555377777189044859520011376295641074225820010046828103529085821048767532072171069130907617986568145052339630099417371229841775083223948249837493631484122886378546453666628015275684870930886802099268366544219399980259811881167537209478988275476808996363584656013453436548981512541129929800863705350208358562474161126930882187788305883374403382713554889961137697
?

[kruoli]

What about something like:

M=7;forprime(a=2,1000,print(a);s=Set([2,a]);for(k=1,M,m=#s;v=vecsort(eval(s));print(k" "#v" "v);print();if(M==k,break);b=1;foreach(v,c,b=b*c);fordiv(b,d,p=d-1;if(ispseudoprime(p),s=setunion(s,Set(p))));if(#s==m,print("No further grow detected, break.");break)))

It will print the first M iterations for every prime up to thousand (if it does not grow, it will break). It runs in less than a second for me.

[Dr Sardonicus]

Nice little routine! Minor nit-pick: Start with a = 3, not 2

I fiddled around with the output commands. I discarded the print() commands, and instead used a vector of size M to record, for each prime a, how many iterations occurred. I had the routine print the primes for which the set hadn't stopped growing when the maximum number of iterations was reached, and after finishing the loop, the vector. Of course, 2*a - 1 is composite for most a.

By extending the range, I got the routine to seemingly stall.

Further investigation showed that for the primes p = 1531 and 6121, the sets grow even faster than for p = 7.

[Dr Sardonicus]

The following code uses an alternative to fordiv(), which seemed to be creating a bottleneck. It aborts if the set gets "too big." It also suppresses output for sets [2, a] when 2*a - 1 is composite. It seems to run fairly quickly.

? c=0;forprime(a=3,2000,s=Set([2,a]);ok=1;k=0;until(ok==0,k++;v=vecsort(eval(s));n=#v;w=vector(#v-1,i,v[i+1]);for(i=1,2^#w - 1, b=binary(i);bpad=vector(#w-#b);b=concat(bpad,b);f=1;for(j=1,#w,if(b[j]==1,f*=w[j]));p=2*f-1;if(ispseudoprime(p),s=setunion(s,Set(p))));if(#s==n&&k==1,c++;break);if(#s==n&&k > 1,v=vecsort(eval(s));print("For p = "a" sequence ends after "k" iterations with "v);ok=0);if(#s > 16,print("For p = "a" terminated after "k" steps with set of size " #s);ok=0)));print("There were "c" odd primes a for which 2*a - 1 was composite.")
For p = 3 terminated after 6 steps with set of size 39
For p = 7 terminated after 6 steps with set of size 29
For p = 19 sequence ends after 3 iterations with [2, 19, 37, 73]
<snip>
For p = 829 sequence ends after 6 iterations with [2, 829, 1657, 3313, 5492953, 10979281, 49995895719789433, 3637137214184962843457821297]
<snip>
For p = 1279 sequence ends after 4 iterations with [2, 1279, 2557, 5113, 26147881]
For p = 1297 sequence ends after 3 iterations with [2, 1297, 2593, 6726241]
For p = 1399 sequence ends after 2 iterations with [2, 1399, 2797]
For p = 1429 sequence ends after 2 iterations with [2, 1429, 2857]
For p = 1459 sequence ends after 2 iterations with [2, 1459, 2917]
For p = 1531 terminated after 5 steps with set of size 34
For p = 1609 sequence ends after 2 iterations with [2, 1609, 3217]
For p = 1627 sequence ends after 2 iterations with [2, 1627, 3253]
For p = 1657 sequence ends after 3 iterations with [2, 1657, 3313, 10979281]
For p = 1759 sequence ends after 2 iterations with [2, 1759, 3517]
For p = 1867 sequence ends after 2 iterations with [2, 1867, 3733]
There were 253 odd primes a for which 2*a - 1 was composite.
?

I initially used a fragment of the above code to see how the set grows for individual primes. The following was spectacular:

? s=Set([2,6121]);v=vecsort(eval(s));n=#v;w=vector(#v-1,i,v[i+1]);for(i=1,2^#w - 1, b=binary(i);bpad=vector(#w-#b);b=concat(bpad,b);f=1;for(j=1,#w,if(b[j]==1,f*=w[j]));p=2*f-1;if(ispseudoprime(p),s=setunion(s,Set(p))));print(#s)
3
? v=vecsort(eval(s));n=#v;w=vector(#v-1,i,v[i+1]);for(i=1,2^#w - 1, b=binary(i);bpad=vector(#w-#b);b=concat(bpad,b);f=1;for(j=1,#w,if(b[j]==1,f*=w[j]));p=2*f-1;if(ispseudoprime(p),s=setunion(s,Set(p))));print(#s)
5
? v=vecsort(eval(s));n=#v;w=vector(#v-1,i,v[i+1]);for(i=1,2^#w - 1, b=binary(i);bpad=vector(#w-#b);b=concat(bpad,b);f=1;for(j=1,#w,if(b[j]==1,f*=w[j]));p=2*f-1;if(ispseudoprime(p),s=setunion(s,Set(p))));print(#s)
8
? v=vecsort(eval(s));n=#v;w=vector(#v-1,i,v[i+1]);for(i=1,2^#w - 1, b=binary(i);bpad=vector(#w-#b);b=concat(bpad,b);f=1;for(j=1,#w,if(b[j]==1,f*=w[j]));p=2*f-1;if(ispseudoprime(p),s=setunion(s,Set(p))));print(#s)
16
? v=vecsort(eval(s));n=#v;w=vector(#v-1,i,v[i+1]);for(i=1,2^#w - 1, b=binary(i);bpad=vector(#w-#b);b=concat(bpad,b);f=1;for(j=1,#w,if(b[j]==1,f*=w[j]));p=2*f-1;if(ispseudoprime(p),s=setunion(s,Set(p))));print(#s)
464
?