Possible values for k (in Mersenne factors)

[science_man_88]

I realized there is something about the k values that work for double mersennes that are different since the single mersenne has form q=2*j*p+1; 2*k*q+1 = 2*k*(2*j*p+1)+1 , which is 2*k+1 mod p so if k is congruent to i mod p where p=2*i+1 the whole thing divides by p and isn't prime. but it's it's probably not all that helpful in codes.

[ATH]

Factors of 2^{(2^p-1)} - 1 have the form q = 2*k*(2^p-1) + 1
and not 2*k*(2*j*p+1) + 1 ?

[science_man_88]

Originally Posted by ATH

Factors of 2^{(2^p-1)} - 1 have the form q = 2*k*(2^p-1) + 1
and not 2*k*(2*j*p+1) + 1 ?

the original mersenne they are based off has form q=2*j*p+1 so 2*k*q+1 expands into 2*k*(2*j*p+1)+1

[Batalov]

Originally Posted by science_man_88

the original mersenne they are based off has form q=2*j*p+1 so 2*k*q+1 expands into 2*k*(2*j*p+1)+1

Yes! It does!
It expands into 2*k*(2*j*p+1)+1, where j has the only one value: j = (2^p-1-1)/(2*p).

And then it contracts to 2*k*(2^p-1)+1.

[science_man_88]

Originally Posted by Batalov

Yes! It does!
It expands into 2*k*(2*j*p+1)+1, where j has the only one value: j = (2^p-1-1)/(2*p).

And then it contracts to 2*k*(2^p-1)+1.

right but in the other form we can show it as 2k+1 mod p which you can reduce to if k is i mod p where p=2*i+1 then the candidate divides by p. like I said it probably won't speed up code ( it's already likely known but it's something the k in the single mersennes don't work off no candidate is 0 mod p there.

[ATH]

Your version 2*k*(2*j*p+1) + 1 is not correct:
https://primes.utm.edu/notes/proofs/MerDiv.html

The factor of a Mersenne number 2^p-1 is of the form 2*k*p + 1.

You are saying the factor is 2*k*q + 1 where q is itself of the form of a mersenne factor of the exponent. But the exponent is prime so the only factor is the number itself, which is why your j variable has only 1 possible value.

[science_man_88]

Originally Posted by ATH

Your version 2*k*(2*j*p+1) + 1 is not correct:
https://primes.utm.edu/notes/proofs/MerDiv.html

The factor of a Mersenne number 2^p-1 is of the form 2*k*p + 1.

You are saying the factor is 2*k*q + 1 where q is not the exponent of the double mersenne prime, but instead a factor of the mersenne prime in the exponent, but this is not correct.

it is !! in MMp Mp is the exponent and Mp takes the form q=2*j*p+1 MMp is of form 2*k*q+1 where q=2*j*p+1 so it's equal to 2*k*(2*j*p+1)+1 = 4*k*j*p+2*k+1 .

[hansl]

Originally Posted by LaurV

Sorry that I am waking up this topic. I searched for it in the light of this new discussion, I did not get the right thing I was searching for, but reaching this topic remembers me that I had implemented the "420 classes" not long after this discussion, doing a pari script to generate the vectors. The code might be useful to somebody, so I will include here. It is about 15% faster than the "60 classes" code described in the posts above. It also eliminates the restriction that the starting q be a multiple of 420. The version with 4620 classes would be slower because of overhead (hey! pari is not mfaktc!).

Somebody can use it for understanding how mfaktc works, or for TF-in small exponents which are outside of the mfaktc range (again, there are much faster programs there outside, this is not for "production" purposes).

Attachment 11028
[can't insert code, is over 10k characters limit]

I know this post is from like 5yrs ago, but I've found this useful and made some improvements to this PARI script that I figured I would share.

First I just generate a vector of vectors for the "classes" and just index to the appropriate one as needed instead of having hundreds or thousands of "if" statements.
Also I have changed the format of how it stores these classes slightly, so it keep the absolute k offset rather than the difference between adjacent ones. This makes it simpler in the main loop to just go through a constant step width and processing one class at a time, which I think makes things a little faster. With this code it is trivial to do 4620 classes, which also gives its own speedup over the previous 420. Even 4620*13 is possible if you raise parisize or parisizemax to about 8GB.

The biggest speedup though is making the computation parallel/multithreaded using parfor. (Note to any Windows users: the default build doesn't seem to support multithreading, I found a special multithreaded build for windows, but it seemed not very effective and I actually got better speeds running PARI under windows subsystem for linux "WSL", with Debian).

\\ "absolute" classes, not delta-based
createclass_abs(n,p,l) = {
if(gcd(p,n)==1,
    my(k=0,z,q);
    vector(l,dummy,
        k++;
        while((gcd(q=2*k*p+1,n)!=1 || (z=q%8)==3 || z==5)&&k<=n,k++); k
    ),
    []
)
}

\\ create all classes, n=4620 recommended in most cases, next level would be 4620*13, which takes ~8G RAM and quite a while to generate
createclasses_abs(n)={ my(l=eulerphi(n)); vector(n-1,p,createclass_abs(n,p,l)) }

\\ main trial factoring function called by par_tf
tf_cl_s(p,fromq=0,toq=2^32,cl,ncl,v,logfile) = {
    my(fromk,q,z,qq,kn,fm,milestone=max(1,ncl\1000),goal=milestone);
    if(#v==0,error("p should better be a prime"));
    kn=max(0,fromq-1)\(2*p*cl)*cl;
    for(i=1,#v,
        fromk=kn+v[i];
        forstep(q=fromk*p<<1+1, toq, cl*p<<1,
            if(Mod(2,q)^p==1,
                printf("M%d has a factor: %d                                                                                 \n",p,q);
                if(logfile,write(logfile,"M"p" has a factor: "q));
            );
        );
        \\if(i>=goal,printf("%3.2f%% done. class %d/%d (%d) fromk=%d %c",100*i/ncl,i,ncl,v[i],fromk,13); goal+=milestone);
    );
}

\\ TF a vector of Mersenne prime exponent candidates, checking for factors in the range: [2^lo_bits,2^hi_bits]
\\ multithreaded TF with logging, with each exponent in a separate thread
par_tf_pv(pv,lo_bits,hi_bits,cl,ncl,vv,logfile) = {
    my(func=tf_cl_s,fromq=2^lo_bits,toq=2^hi_bits); 
    print("processing ",#pv," primes in list");
\\ "vv[pv[i]%cl]" gives use the "classes" vector we want without hundreds of if statements
    parfor(i=1,#pv,func(pv[i],fromq,toq,cl,ncl,vv[pv[i]%cl],logfile));
}

\\ Example TF of a range of exponents:
\\ cl = 4620;
\\ vv = createclasses_abs(cl);
\\ par_tf_pv(primes([10^9,10^9+10^6]),1,55,cl,eulerphi(cl),vv,"test.log")


\\ For multithreaded TF of a single exponent, with each class in its own thread
par_tf_cl_s(p,fromq=0,toq=2^32,cl=420,ncl,v,logfile) = {
    my(fromk,q,z,qq,kn,fm,milestone=max(1,ncl\10000),goal=milestone);
    if(#v==0,error("p should better be a prime"));
    kn=max(0,fromq-1)\(2*p*cl)*cl;
    parfor(i=1,#v,
        fromk=kn+v[i];
        forstep(q=fromk*p<<1+1, toq, cl*p<<1,
            if(Mod(2,q)^p==1,
                printf("M%d has a factor: %d                                                              \n",p,q);
                if(logfile,write(logfile,"M"p" has a factor: "q));
            );
        );
        if(i>=goal,printf("%3.2f%% done. class %d/%d (%d) fromk=%d %c",100*i/ncl,i,ncl,v[i],fromk,13); goal+=milestone);
    );
}

\\ Example TF of single exponent in parallel (note this version takes literal low/high factor limits, not bit levels)
\\ cl = 4620*13*17;
\\ ncl = eulerphi(cl);
\\ e = 1249;
\\ par_tf_cl_s(e,2^1,2^32,cl,ncl,createclass_abs(cl,e%cl,ncl));

I created two main functions for parallel TFing. The first parallelizes across a vector of exponents(for testing many exponents to relatively low bit depth), and the other is for TF of a single exponent, which parallelizes the computation across each class.

I added some comments which hopefully make sense, and a simple example of how each version can be used. When doing a single exponent you can actually go much higher in terms of # of classes computed, since it only needs to compute a single set of classes for that exponent.

Oh and lastly my version of the code doesn't stop on the first factor found, but does the whole range anyways. This should be trivial to change by adding a return or break if you want though.

Edit: I forgot to mention that you can also set "nbthreads" default in "~/.gprc" config file, or using the command: "default(nbthreads, x)" in gp prompt, to determine how many threads are used by parfor. It seems to default to all cores + hyperthreads.

[LaurV]

Nice to see that people still play with this stuff. We salute the parallelizing (ha! this is actually a word, it is not underlined with red!) of the code, as well as the arbitrary class extension. Unfortunately we are at job now (lunch break) and can't really play, but we are eager to go home tonight and give it a spin. We wanted to make a parallel TF version but never found the time. As pari does not care much about the size of the exponents, this may be quite useful for James&co to extend the 10G factorization, if it proves reasonable fast. Mind that for mfaktX the exponent size is still limited to 32 bits.

[hansl]

I'm tempted to implement this code directly in C++/GMP, although that alone probably won't make much improvement since PARI is already running on GMP.
I could try adding some proper sieving on top which would surely help for deeper bits.

Although at that point maybe I'm just re-implementing factor5. I should probably look a little closer at that code...

[hansl]

Originally Posted by LaurV

Unfortunately we are at job now (lunch break) and can't really play, but we are eager to go home tonight and give it a spin.

Hi LaurV, just curious if you ever got around to playing with the code and if you had any other thoughts?

Originally Posted by LaurV

As pari does not care much about the size of the exponents, this may be quite useful for James&co to extend the 10G factorization, if it proves reasonable fast.

Yep, you called it! This is what I was using verify/complete the factorization up to 55 bits for James' database. I was being a bit coy about it originally before the work was complete, so as not to overwhelm James' servers with results, if others had started submitting the same data etc. (There was a small part of the script that I omitted which was automatically submitting 1M exponent ranges at a time to James' bulk results form)

Now that I'm done with those smaller/faster bit ranges, I think its much less concern, so I'm looking at possibly improving efficiency once again, and seeing how much further we can go with bit-levels. If others want to join that would be nice too! I can start a separate thread with more detail if there is any interest in my crazy sub-project.

I'm doing 55-56 bits right now, sort of working downward from 8G and 9G ranges. I'm leaving more computers and cores assigned to "regular" primenet work than my initial run, so its even slower going this time.

There was one thing I had noticed from your earlier attachment: https://www.mersenneforum.org/attach...chmentid=11028
If we look at your matrix of classes mod 420, you can see that each row is actually the same as all the others, just rotated by some amount (its all basically a single wheel).
For example if you select a subset of say 1/3rd of the values from any given row in a text editor and search/highlight other occurrences, you can see most of the other rows get highlighted (really all of them match, but some wrap around from end to beginning, so text editor obviously wouldn't catch that). You can even select all but the first or last value from any given row, and there will be one other corresponding row shifted by one.

So right now I hypothesize that there exists some simple(without any loops) modular arithmetic to determine the amount of rotation/offset between any two rows, given their residues, but I haven't worked it out yet. That would make it much simpler to only store one row(the first one, for example) and offset into it by the appropriate amount.

Any ideas?