20180622, 20:59  #1 
"Robert Gerbicz"
Oct 2005
Hungary
1,409 Posts 
Saving tons of PRPC and (hypotetical) LLC tests, a new method
PRPC refers to (weak) Fermat test for (2^p1)/d.
And LLC refers to a LL type test for (2^p1)/d, see http://mersenneforum.org/showthread....732#post471732 ofcourse currently the main issue with LL/LLC is that we have no strong error check for these LL (type) tests. But first concentrating on PRPC tests. You can see such lines when you can check out the recent results page: https://www.mersenne.org/report_expo...1344943&full=1 https://www.mersenne.org/report_expo...2389403&full=1 https://www.mersenne.org/report_expo...2462281&full=1 https://www.mersenne.org/report_expo...3009703&full=1 https://www.mersenne.org/report_expo...4373911&full=1 What happens here? We do a PRPC test on (2^p1)/d1, but when a new prime divisor appears, then we redo this, because the full residue is not available. Ofcourse there is a folklore solution: store the full residue. But that takes lots of storage and bandwidth. I can do much better here, and this is the trick. Let p>3 prime and mp=2^p1. Compute Code:
[1]: 3^mp==res mod mp, independently from d. If mp/d is prime, then using Fermat's theorem we know that: 3^(mp/d)==3 mod (mp/d), so Code:
[2]: 3^mp==3^d mod (mp/d) also holds (note that here we're already weaker than Fermat, but not much weaker). 3^mp==res==3^d mod (mp/d) let s=3^d mod (mp/d), then Code:
[3]: res=s+w*(mp/d), for some w integer, where 0<=w<d. Multiply by d, we get d*(ress)=w*mp, see this equation mod 2^t, where t<p, then w==d*(sres) mod 2^t, but if 2^t>d, then we have a unique candidate for w: Code:
[4]: w=(d*(sres) mod 2^t) (not forget the only conditions: t<p and d<2^t) It was a quite general idea, works for all numbers, just go back to [3], suppose N is odd, and d is a nontrivial divisor of N. [3]': res=s+w*(N/d), see this mod 2^t w==(ress)*inv mod 2^t, where (N/d)*inv==1 mod 2^t, here inv exists, because N/d is also odd. and if d<2^t and d<=w then N/d is composite. [ where w=((ress)*inv) mod 2^t) ]. Not get confused, but in Mersenne number case it was inv=d for p>t ! A little quick and dirty test, just get the small q prime divisors for mp where q<10000*p , and test mp/d, where d is the product of all found prime divisors of mp. Code:
prp(n)=return(Mod(3,n)^n==Mod(3,n)) prpc(maxp)={t=512;num_large_d=0;forprime(p=t,maxp, d=1;mp=2^p1;forstep(q=2*p+1,10000*p,2*p,if(isprime(q)&&lift(Mod(2,q)^p)==1,d*=q)); if(d>1&&d<mp,res=lift((Mod(3,mp)^mp));s=lift(Mod(3,mp)^d)%(mp/d); w=(d*(sres))%(2^t); if(d>2^t,num_large_d+=1); if(prp(mp/d),print(["Found prp ",p,d])); if(w<d&&2^t>d,print(["Candidate! ",p,d])))); print("Number of large divisors=",num_large_d);} prpc(10000) only s has a slightly complicated formula. Code:
prp(n)=return(Mod(3,n)^n==Mod(3,n)) fmult(u,v,p)=return([(u[1]*v[1]+3*u[2]*v[2])%(2^p1),(u[1]*v[2]+u[2]*v[1])%(2^p1)]) fpowmod(v,e,p)={ ret=[1,0];h=v; while(e, if(e%2,ret=fmult(ret,h,p)); h=fmult(h,h,p);e>>=1); return(ret)} LLT(n,x0,it)={mp=2^n1;x=Mod(x0,mp);for(i=1,it,x=x^22);return(lift(x))} llc(maxp)={t=512;forprime(p=t,maxp, d=1;mp=2^p1;forstep(q=2*p+1,10000*p,2*p,if(isprime(q)&&lift(Mod(2,q)^p)==1,d*=q)); if(d>1&&d<mp, res=LLT(p,4,p); s=(2*(fpowmod([2,1],d+kronecker(3,mp/d),p)[1]))%(mp/d); w=(d*(sres))%(2^t); if(d>2^t,num_large_d+=1); if(prp(mp/d),print(["Found prp ",p,d])); if(w<d&&2^t>d,print(["Candidate! ",p,d])))); print("Number of large divisors=",num_large_d);} llc(10000) Code:
["Found prp ", 881, 26431] ["Candidate! ", 881, 26431] ["Found prp ", 1459, 1362463807] ["Candidate! ", 1459, 1362463807] ["Found prp ", 3041, 135391639199] ["Candidate! ", 3041, 135391639199] ["Found prp ", 3359, 6719] ["Candidate! ", 3359, 6719] ["Found prp ", 4127, 154517628583] ["Candidate! ", 4127, 154517628583] ["Found prp ", 4243, 101833] ["Candidate! ", 4243, 101833] ["Found prp ", 7673, 2563193011919] ["Candidate! ", 7673, 2563193011919] Number of large divisors=0 ps. say if 2^t>d>2^(t1) then you'll have a good chance that w<d, but this doesn't mean that you have found a prime, because with more than 50% chance w<d will be true. Note that for almost all known smallish p value (p<10^9 or so) it is true that d<2^512. So you'd only need additional PRPC tests for those d>2^512 and for those that survived this new test. For this we'd need to store RES512, if we'd choose t=512. And not forget that with this mod 2^t shrink we are (still) weaker than a standard Fermat test, but not much weaker as you can see from the example runs. 
20180623, 16:30  #2 
"TF79LL86GIMPS96gpu17"
Mar 2017
US midwest
2^{2}×3^{2}×127 Posts 
Great! And the questions begin
Excellent. While this does not seem to impact the primalityunknown PRP or LL, as I incompletely understand it, it should really help the PRPC and LLC progress. Like your previous insights, it will take a while to understand, implement, roll out, and yield advantages. And that long process starts when you introduce such things, so keep them coming!
So now some questions, A) step one is the PRP test if no factors are found yet. Always save res512 from that? B) does the first cofactor run generate the res512 and send it to the (enhanced to accept it) PrimeNet server? C) what's the language of your code sections for PRPc and llc? D) For confidence in accurate transmission of res512 values, do we use some agreed standard like CRC32 from the start? E) Is addition of res512 storage practical on the server? A question for Madpoo et al F) What's inv there? Inverse element in some group? "Not get confused" (Too late...;) Last fiddled with by kriesel on 20180623 at 16:33 
20180623, 19:22  #3  
"Robert Gerbicz"
Oct 2005
Hungary
1,409 Posts 
Quote:
Quote:
Quote:
number of factors of mp. And yes, we'd need to store it somewhere. Quote:
Quote:
It's more than a group, it is a ring, but yes ofcourse inv is the multiplicative inverse of N/d in Z[2^t]. In the general case there is (ress) factor in the "final" formula, and in the Mersenne case there was sres in [4], a negative sign came from inv, because inv=d. Returning to LL case: since we need p iterations here and we can't get even the low 64 bits of S[p] if you know (S[p2] mod 2^64), we really need that two (cheap) LL iterations, and we need to keep those S[p2] residues also, because: 1. for double checking the old tests 2. S[p]=2 is necessity if mp is prime, but I don't see the sufficiency [it could be true, but this needs a proof]. 

20180623, 19:31  #4  
"Robert Gerbicz"
Oct 2005
Hungary
1,409 Posts 
Quote:
Code:
[4]: w=(d*(sres') mod 2^t) Not written, but a general trick works here if you used 3^(mp+1) or 3^(mp1) instead of 3^mp, then there are 3 possibilities for res'. Just for completeness, the two corrected tests, when we store the shrinked residues: Code:
prp(n)=return(Mod(3,n)^n==Mod(3,n)) prpc(maxp)={t=512;num_large_d=0;forprime(p=t,maxp, d=1;mp=2^p1;forstep(q=2*p+1,10000*p,2*p,if(isprime(q)&&lift(Mod(2,q)^p)==1,d*=q)); if(d>1&&d<mp,res=(lift((Mod(3,mp)^mp)))%(2^t);s=lift(Mod(3,mp)^d)%(mp/d); w=(d*(sres))%(2^t); if(d>2^t,num_large_d+=1); if(prp(mp/d),print(["Found prp ",p,d])); if(w<d&&2^t>d,print(["Candidate! ",p,d])))); print("Number of large divisors=",num_large_d);} prpc(10000) Code:
prp(n)=return(Mod(3,n)^n==Mod(3,n)) fmult(u,v,p)=return([(u[1]*v[1]+3*u[2]*v[2])%(2^p1),(u[1]*v[2]+u[2]*v[1])%(2^p1)]) fpowmod(v,e,p)={ ret=[1,0];h=v; while(e, if(e%2,ret=fmult(ret,h,p)); h=fmult(h,h,p);e>>=1); return(ret)} LLT(n,x0,it)={mp=2^n1;x=Mod(x0,mp);for(i=1,it,x=x^22);return(lift(x))} llc(maxp)={t=512;forprime(p=t,maxp, d=1;mp=2^p1;forstep(q=2*p+1,10000*p,2*p,if(isprime(q)&&lift(Mod(2,q)^p)==1,d*=q)); if(d>1&&d<mp, res=LLT(p,4,p)%(2^t); s=(2*(fpowmod([2,1],d+kronecker(3,mp/d),p)[1]))%(mp/d); w=(d*(sres))%(2^t); if(d>2^t,num_large_d+=1); if(prp(mp/d),print(["Found prp ",p,d])); if(w<d&&2^t>d,print(["Candidate! ",p,d])))); print("Number of large divisors=",num_large_d);} llc(10000) Last fiddled with by R. Gerbicz on 20180623 at 19:38 Reason: adds+grammar 

20180623, 19:48  #5 
"Forget I exist"
Jul 2009
Dumbassville
20261_{8} Posts 
it needs to fall back to the 2 case eventually because if not you won't get anything but a list of 4's. The only other possible way around this is the previous mersenne's square (congruent to the power of 2 next door mod the higher mersenne) or another pair of a and a having their squares equal to 4 mod once you fall to the 2 path backwards you hit 0 mod as a square and sqrt.
Last fiddled with by science_man_88 on 20180623 at 19:49 
20180624, 14:43  #6  
"Robert Gerbicz"
Oct 2005
Hungary
2601_{8} Posts 
Quote:
Ofcourse even with a proof as stated above for obvious reason(s) we'd need to store the RES64 for S[p2]. For my algorithm this question isn't that very interesting (at database or in other level), just noted as a somewhat nice open(?) question. 

20180627, 05:43  #7 
"Mihai Preda"
Apr 2015
2^{2}·3^{2}·5·7 Posts 
Would it be useful/desirable for LL and PRP software to report to the server a larger 512bit residue instead of the current 64bit?

20180627, 08:29  #8  
Jun 2003
2^{2}·11·107 Posts 
Quote:
But if you're going this route, might as well futureproof it and go with a bigger RES1024 or something ( this will fit in a varbinary(128) in sql server, so storage cost isn't high ). 

20180627, 09:22  #9  
"Mihai Preda"
Apr 2015
2^{2}×3^{2}×5×7 Posts 
Quote:
Also, I can do either 512 or 1024. I'd like to know what George thinks. If he gives the green light I can do that. (also we'd need to chose a key for reporting it in the JSON. Probably in addition to the current res64, in a first phase). Last fiddled with by preda on 20180627 at 09:22 

20180627, 09:28  #10 
"Mihai Preda"
Apr 2015
2^{2}·3^{2}·5·7 Posts 
It would be nice if we had a way to derive the current res64 from the new res512 or res1024. My feeling is that it could be done with some guard highbits ("before 0") that would sink the carry propagation from the div3 or mul3.
Or, could the new "bigresidue" be not starting at bit0, but e.g. centered on bit0? i.e. shifted some amount such that it contains bits on both sides of the current res64. Last fiddled with by preda on 20180627 at 09:35 
20180627, 10:26  #11  
"Robert Gerbicz"
Oct 2005
Hungary
10110000001_{2} Posts 
Quote:
We could switch to res1024 only for say p>2^32. And ofcourse you'd to store it only if you have a completed LL/PRP test to save space in database. Quote:
Code:
let f(p,e)=mp=2^p1;return(lift(Mod(3,mp)^e)%(2^512)) then f(p,e+1)=(3*f(p,e)c*mp)%(2^512), but mp==1 mod 2^512, for p>512. so we can write: Code:
f(p,e+1)=(3*f(p,e)+c)%(2^512) for some c=0,1,2 and if you solve it for f(p,e) then: f(p,e)=lift(Mod((f(p,e+1)c),2^512)/3) for some c=0,1,2 (after the additional mod 2^64), ofcourse you would have 3 residues even for mod 2^512. The situation is similar if you have general number (not Mersenne), or base!=3 or p<512. (but that is really boring, since those are factored), or you have f(p,e+2) and want f(p,e); but in this case you can also just iterate the above method, you get 9 possible residues for f(p,e). And again you can say in the case of LL test nothing from RES512 at iteration p, due to the nonlinear squaring. [That is why we'd need to store also the "standard" RES64 at iteration (p2)]. Last fiddled with by R. Gerbicz on 20180627 at 10:36 

Thread Tools  
Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Daylight Saving Time  davieddy  Lounge  27  20150226 18:45 
Prime95 NOT Saving!!!  Primeinator  Information & Answers  9  20080922 19:00 
Not Saving  BioRules  Information & Answers  9  20080531 13:52 
Saving computation in ECM  dave_dm  Factoring  8  20040612 14:18 
saving over a network  crash893  Software  11  20040506 14:15 