Sieving freakishly big MMs (was "World record" phone number?)

[ewmayer]

Originally Posted by ET_

Hi Ernst, the last 4 DM have been sieved up to 7G.

Should we proceed with PARI up to a defined level, or you may show us a better sieving method?

Luigi

I hacked together a deeper version of my factor.c sieve (which is limited to primes < 2^32) into a small standalone code yesterday, and am doing some deeper sieving of the small-k candidates. Sample timing: to sieve q = 2.k.M(43112609) +1 to depth 10^10 needs 26 min on a single 2GHz core of my vintage-2009 macbook. How does that compare to what other folks are using?

(I reproduced the small factor for k = 5; k= 185 I have sieved to 10^11, no factors below that bound).

I am working on some optimizations today which I hope will double (or more) the speed; once I'm satisfied it's as fast as I can make it without major (week or more) effort I'll post it here.

[science_man_88]

Originally Posted by LaurV

You have (a) q=2kp+1 and either (b) q=8x+1 or (c) q=8x-1. Combining (a) and (b) gives q=2*(4z)*p+1 (i.e. k=0 mod 4) in all cases, regardless of p. Combining (a) and (c) can't be solved in general case, you have to split the two cases where p=1 (mod 4) or p=3 (mod 4). First gives k=3 (mod 4) and second gives k=1 (mod 4).

Therefore, factors of M(4z+3) are either 8kp+1 or 8kp+2p+1 (that is 2*(4x)*(4z+3)+1 or 2*(4x+1)*(4z+3)+1 with natural x).
Ex: M11 has factors 23 and 89, with k=1=4*0+1, and k=4=4*1 (renaming the right k). M43 has factors 431, 9719 and 2099863, with correspondent k: 5=4*1+1, 113=4*28+1, 24417=4*6104+1. The k's for such numbers can only be 1,4,5,8,9,12,13, etc.

Same, factors of M(4z+1) are either 8kp+1, or 8kp+6p+1 (that is 2*(4x)*(4z+1)+1 or 2*(4x+3)*(4z+1)+1 with natural x. Ex: M29 has factors 233, 1103 and 2089, with correspondent k's 4=4*1, 19=4*4+3, 36=4*9. M37 has factor 223 with k=3=4*0+3.

Other combinations give composite numbers or primes which are 3 or 5 mod 8, as Batalov said. Those can't be "prime factors" of Mp with prime p.

As Mp is always 3 mod 4, it follows that if Mp is prime, MMp may have only factors with k=1 or k=0 mod 4.

Remark also that the "plus" factors (the one of the form 1 mod 8) can be in any number, including zero, but the "minus" factors (the one being "-1" mod 8) must be in odd number, because Mp is 3 mod 4, and product of two numbers which are 3 (mod 4) is 1 (mod 4). So it can't be an "even" number of "minus" factors. This implies that if Mp is composite, then it has al least one proper factor of the "minus" type (whose k is either 1 or 3 mod 4, depending on p).

To filter (I can't say sieve, but very elementary) possible factors for MMp you can try this small pari function to get a better insight.

mmp47(k,startq,stopq)=
{
    p=43112609;
    until(q>stopq,
        if(k%4,k+=3,k++);
        cnt=0;
        q=startq-1;
        while(q<=stopq && 2*k*(Mod(2,q=nextprime(q++))^p-1)+1,
            if(q>cnt,printf("... %d ...%c",q,13);
            cnt=cnt+10^6)
        );
        print(k" : "q"                 ");
    );
    return(q<=stopq);
}

call it with "mmp47(0,0,10^10)" and wait... you can stop it with ctrl+c when it reaches 185, and restart it again with "mmp47(185,0,10^10)" to skip to the next value after 185 (which would be 188). Next "hump" is 201.

@sm88: I got your PM, I hope the explanations above satisfy you too. Also, 1G4=1.4*10^9 (from electronics, resistors with 1.27 kilo-ohms are marked 1k27, to save space on small PCBs, components, etc, generally the notation is used in physics, 3M5 means 3.5 mega, or 3500000, etc, like rounding where the "small amounts" are not interesting). The k=185 is tested to much higher values, you can use the code above but is slow, some guys tested it to 1T or so (I won't bet my life on it, take it like rumors only).

(17:34)>for(x=1,43112609,if(43112609%x==(x-1),print(x)))
1
2
3
5
6
9
10
15
18
30
45
90
479029
958058
1437087
2395145
2874174
4311261
4790290
7185435
8622522
14370870
21556305

I'm guessing all k values 1 mod (2^(this list)-1) have been taken off ? I say that because (2^43112609-1) mod (2^x-1)=(2^(x-1)-1) and k=1 mod (2^x-1) gives 2*1*(2^(x-1)-1)+1 = (2^x-1) mod (2^(x-1)-1)

edit: the last mod is listed wrong it's mod (2^x-1)

[ET_]

Originally Posted by ewmayer

I hacked together a deeper version of my factor.c sieve (which is limited to primes < 2^32) into a small standalone code yesterday, and am doing some deeper sieving of the small-k candidates. Sample timing: to sieve q = 2.k.M(43112609) +1 to depth 10^10 needs 26 min on a single 2GHz core of my vintage-2009 macbook. How does that compare to what other folks are using?

(I reproduced the small factor for k = 5; k= 185 I have sieved to 10^11, no factors below that bound).

I am working on some optimizations today which I hope will double (or more) the speed; once I'm satisfied it's as fast as I can make it without major (week or more) effort I'll post it here.

Thank you very much Ernst! I stopped searching k=185 at 5*10^10...
Your timings seems at least 10x faster than our PARI script: when the code will be ready I think we may arrange a distributed sieving on all k.

I also tested a primality testing of 2*185*M47+1 with pfgw. It is slow, but can be done (it takes about the time of a LL test on a 56M exponent).

I was doing some tests, and I know that LaurV is updating a table for all k candidates below 1000 for MM from 34 to 47.

So, if you all agree, we may wait for the table, define a criterion for the reservations, and then "Gentlemen, start our sievers!" (that should be "your siever", but this way sounds better)...

Luigi

[axn]

Originally Posted by ewmayer

Sample timing: to sieve q = 2.k.M(43112609) +1 to depth 10^10 needs 26 min on a single 2GHz core of my vintage-2009 macbook.

Assuming you have a Core 2, this is at least an order-of-magnitude slower than the state-of-the-art.

NewPGen can sieve a fixed-k range at about 6G / min on similar machine (while it is not exactly apples-to-apples, the complexity of calculation is comparable).

My sieve in PARI can do a 4G range in appr 2 mins. [k=1..10000, primelimit set to 4G]. Since I rely on precomputed primelist, that is as high as I can go with 32-bit PARI. Nonetheless, the point is that you're way off in terms of what is achievable. Could there be an algorithmic issue with your code?

[LaurV]

@sm88: You are a bit confuse here, there is no shame, I was doing this mistake in the beginning, and I am still doing it from time to time... Generally when you switch to modular things, the affirmation "if a=b then x^a=x^b" is not true. There is something called "eulerphi" in the middle... For example, 5=16 (mod 11), but 2^5=32=10 (mod 11) and 2^16=65536=9 (mod 11). When you do double mersenne, i.e MMp, the Mp is the one needing to have those modular properties, and not p. It does not help us too much what properties p (mod x) has; only what properties Mp (mod x) has, helps. You can't "raise" or "lift" from p to Mp. Those are tremendous numbers then, difficult to work with.

[LaurV]

Originally Posted by ET_

I was doing some tests, and I know that LaurV is updating a table for all k candidates below 1000 for MM from 34 to 47.

Yes, and this is what I came with, after a day and a night, with a laptop working unsupervised at home during I was scratching my head and arguing with colleagues at job (about different things).

mmpp.zip

Grrr, I tried to make some snaps, but either the file is too big, the picture is too wide, I have too many attachments (is there no way to overcome that restriction??? xyzzy? I'm not the kind of guy who attach thousands of things to posts just for occupying the space, and generally, beside of my long-speech posts, I don't waste space here!) so I made a zip. Is not nice looking, and you have to unzip it by yourself.

I put the excel table and the pari scripts that helped me, updated. There was a laptop which helped me too, but that I can't attach it to the post (I tried! )

Remark that there are some hidden rows, you can unhide them to see more information, but nothing new, everything is from the web page or wikis, and also remark that I am working the dark-green (khaki?) rows, in fact #39 and #41 should be finished to 20G by now, I let two workers in the morning running mmppsk() for the #k39 and #k41 lists, as shown in the sample code (the laptop is a two-core machine).

Also, most of the work was dupli-tripli-cated, as some people also did the testing, I didn't know ET is interesting in playing with those k's and he didn't know I am spending my time into it, so before we got contact on PM, we did the same work.

edit: some supermod can move this thread at all to the "doppi" subforum (!?) thanks
edit2: obviously, I use the notation M#34 to M#47, like "number", "index", "ordinal", "the 34th or 47th mersenne prime number", to make distinction between M#6 (which would be the sixth mersenne prime number, M#xx is always a prime number) and M6 which is 2^6-1. In this light, M31 is 2147483647, but M#31 is 2^216091-1, a prime, and MM31 is a 646456993-digits number. Now someone tells me how many digits in MM#31, and how many in M#M31 and respective in M#M#31, and if they are prime...

[ET_]

Originally Posted by LaurV

Yes, and this is what I came with, after a day and a night, with a laptop working unsupervised at home during I was scratching my head and arguing with colleagues at job (about different things).

Attachment 8674

Grrr, I tried to make some snaps, but either the file is too big, the picture is too wide, I have too many attachments (is there no way to overcome that restriction??? xyzzy? I'm not the kind of guy who attach thousands of things to posts just for occupying the space, and generally, beside of my long-speech posts, I don't waste space here!) so I made a zip. Is not nice looking, and you have to unzip it by yourself.

I put the excel table and the pari scripts that helped me, updated. There was a laptop which helped me too, but that I can't attach it to the post (I tried! )

Remark that there are some hidden rows, you can unhide them to see more information, but nothing new, everything is from the web page or wikis, and also remark that I am working the dark-green (khaki?) rows, in fact #39 and #41 should be finished to 20G by now, I let two workers in the morning running mmppsk() for the #k39 and #k41 lists, as shown in the sample code (the laptop is a two-core machine).

Also, most of the work was dupli-tripli-cated, as some people also did the testing, I didn't know ET is interesting in playing with those k's and he didn't know I am spending my time into it, so before we got contact on PM, we did the same work.

edit: some supermod can move this thread at all to the "doppi" subforum (!?) thanks

Thank you LaurV.

1 - I will update the history page and (try to) arrange a new subpage of the site related to sieving (but not before next week).

2 - Meanwhile, if no one else will, I am working on MM47, leaving MM40 and MM44 to volunteers for further sieving.

3 - Ernst and axn may test their sievers against LaurV script, and let us know about correlated efficiency.

4 - My idea for the next week(s) would be to deliver some web-service to:
authenticate a user
contact the db
download a k and its limit
upload k with the the new limit
check the result's coherence and update the db
Between the upload and the download there should be some program/script apt to sieve, and maybe releasing an output verifiable by the server.
We may set limits to the sieving job, and once all candidates have been taken to that limit, decide whether raise the limit, start LLtesting or whatever...

Does the idea interest anyone?

Luigi

[axn]

Originally Posted by ET_

Does the idea interest anyone?

If you're gonna do it, you might as well do it right.

1. Sieving is much more efficient than "one k at a time" approach. In fact, if we can find an optimal (or pretty good) addition/subtraction chain for all the exponents in question, sieving all the MMps together will be the most efficient. However, the sieving approach would need us to specify the k-range of interest (k < 10^6?), such that P95 can efficiently support them (i.e. have a good FFT to test them).

2. None of the scripts/program so far are even close to optimal. If you're gonna throw compute power at it, perhaps it'll pay to invest in some upfront programming. I know that Ken_g6 has tackled k*2^n+1 type sieving -- perhaps he can adapt the code for our purpose? [I would anticipate an optimal siever to do a range of 10^12/day/core -- so no real advantage in pushing thru with the current scripts]

3. Once we sieve something to a reasonable depth (10^15?), we can then P-1 and/or ECM them before doing PRP. Since these would be very large primes on their own right, PRP before MMp trial division sounds right.

[science_man_88]

Originally Posted by LaurV

@sm88: You are a bit confuse here, there is no shame, I was doing this mistake in the beginning, and I am still doing it from time to time... Generally when you switch to modular things, the affirmation "if a=b then x^a=x^b" is not true. There is something called "eulerphi" in the middle... For example, 5=16 (mod 11), but 2^5=32=10 (mod 11) and 2^16=65536=9 (mod 11). When you do double mersenne, i.e MMp, the Mp is the one needing to have those modular properties, and not p. It does not help us too much what properties p (mod x) has; only what properties Mp (mod x) has, helps. You can't "raise" or "lift" from p to Mp. Those are tremendous numbers then, difficult to work with.

you clearly forgot the formula for Mp = 2^(p-y)*My+(2^(p-y)-1) this shows the relationship I talked about it shouldn't be a mistake as we are dealing with 2*k*p+1 and p in this case is a mersenne number. under this formula, "if a=b then x^a=x^b" is true.

edit:actually (x^a-1)=(x^b-1) x=2 is true sorry.

[LaurV]

Originally Posted by axn

Sieving is much more efficient than "one k at a time" approach
<snip>
I know that Ken_g6 has tackled k*2^n+1 type sieving -- perhaps he can adapt the code for our purpose?

Of course sieving all would be faster than taking them one by one, but we have no software yet. It would only need a little bit of tickling with NewPGen, that has already k*2^n+1 implemented, it would need only a small addition for k*2^n+1-2k (db)

On the bright side, when taking #39 to 10G, k=225 popped out a factor : 7333051961. I am happy about it, I want to get rid of that long line in the table, to make the table smaller

[LaurV]

Originally Posted by axn

3. Once we sieve something to a reasonable depth (10^15?), we can then P-1 and/or ECM them before doing PRP. Since these would be very large primes on their own right, PRP before MMp trial division sounds right.

This bothers me since I have read it. Especially the last part. "PRP before MMp trial division sounds right". Does it?




(gotcha back )