mersenneforum.org

mersenneforum.org (https://www.mersenneforum.org/index.php)
-   GpuOwl (https://www.mersenneforum.org/forumdisplay.php?f=171)
-   -   gpuOwL: an OpenCL program for Mersenne primality testing (https://www.mersenneforum.org/showthread.php?t=22204)

Prime95 2017-09-21 03:07

So gpuOwl is printing non-standard residues for a Fermat test. That's not particularly good.

Here's my proposal. A Fermat probable prime test for N is 3^(N-1) = 1 mod N. Iteration 0 is 3. Using a standard left-to-right exponentiation algorithm, we get iteration 1 = 3^(top 2 bits of N), iteration 2 = 3^(top 3 bits of N), and so forth.

I believe this is what LLR, PFGW, and prime95 all do for interim and final residues.

You can still do the no-mul-by-3's implementation. That is, N = 2^P-1, and you will calculate 3^(N+1) mod N. To print a standard iteration 1 residue, output 3^(top 3 bits of N+1)/3, etc. That is, you'll do divide-by-3 to output interim residues and divide-by-9 to get the final residue. This is what I'm doing in the next prime95 release -- testing it now.

axn 2017-09-21 03:35

[QUOTE=Prime95;468233]So gpuOwl is printing non-standard residues for a Fermat test. That's not particularly good.

Here's my proposal. A Fermat probable prime test for N is 3^(N-1) = 1 mod N. Iteration 0 is 3. Using a standard left-to-right exponentiation algorithm, we get iteration 1 = 3^(top 2 bits of N), iteration 2 = 3^(top 3 bits of N), and so forth.

I believe this is what LLR, PFGW, and prime95 all do for interim and final residues.

You can still do the no-mul-by-3's implementation. That is, N = 2^P-1, and you will calculate 3^(N+1) mod N. To print a standard iteration 1 residue, output 3^(top 3 bits of N+1)/3, etc. That is, you'll do divide-by-3 to output interim residues and divide-by-9 to get the final residue. This is what I'm doing in the next prime95 release -- testing it now.[/QUOTE]

What all types of test are you planning to cover under the "Gerbicz error check" in P95? Only LL or the general k*b^n+/-1 form? Because, your approach will not work for bases other than 2.

preda 2017-09-21 03:40

[QUOTE=Prime95;468233]So gpuOwl is printing non-standard residues for a Fermat test. That's not particularly good.

Here's my proposal. A Fermat probable prime test for N is 3^(N-1) = 1 mod N. Iteration 0 is 3. Using a standard left-to-right exponentiation algorithm, we get iteration 1 = 3^(top 2 bits of N), iteration 2 = 3^(top 3 bits of N), and so forth.

I believe this is what LLR, PFGW, and prime95 all do for interim and final residues.

You can still do the no-mul-by-3's implementation. That is, N = 2^P-1, and you will calculate 3^(N+1) mod N. To print a standard iteration 1 residue, output 3^(top 3 bits of N+1)/3, etc. That is, you'll do divide-by-3 to output interim residues and divide-by-9 to get the final residue. This is what I'm doing in the next prime95 release -- testing it now.[/QUOTE]

Apparently what I compute is a Strong PRP [url]https://primes.utm.edu/prove/prove2_3.html[/url]
which is stronger than a Fermat PRP. Given no additional cost and some benefit, then why not do the stronger test?

I do not calculate 3^(N+1) mod N, but 3^((N+1)/2) mod N. I think that's exactly the difference that makes it "strong", by halving the probability of false positive.

Prime95 2017-09-21 04:08

[QUOTE=axn;468235]What all types of test are you planning to cover under the "Gerbicz error check" in P95? Only LL or the general k*b^n+/-1 form? Because, your approach will not work for bases other than 2.[/QUOTE]

I'll support all (k*b^n+c)/known_factors.

You are correct, the no-mul-by-3's code is only applies to Mersenne numbers with no known factors. In all other cases, I'll do a standard Fermat PRP test with Gerbicz recovery.

axn 2017-09-21 04:18

[QUOTE=Prime95;468237]I'll support all (k*b^n+c)/known_factors.

You are correct, the no-mul-by-3's code is only applies to Mersenne numbers with no known factors. In all other cases, I'll do a standard Fermat PRP test with Gerbicz recovery.[/QUOTE]

Gerbicz method can only be applied to a sequence of powerings (a^(b^N)), so standard PRP is not compatible with it. Which is why gpuOwl is implementing 3^(N+1) so that the computation sequence is an uninterrupted sequence of squarings only. Now, as luck would have it, we can flip back and forth between standard PRP and powerings-only mode with a simple division/multiplication of 3 (in the case of base 2 only). But that will not work in case of other bases, nor will it work in the case where there is a factor.

preda 2017-09-21 04:33

[QUOTE=Prime95;468233]So gpuOwl is printing non-standard residues for a Fermat test. That's not particularly good.

Here's my proposal. A Fermat probable prime test for N is 3^(N-1) = 1 mod N. Iteration 0 is 3. Using a standard left-to-right exponentiation algorithm, we get iteration 1 = 3^(top 2 bits of N), iteration 2 = 3^(top 3 bits of N), and so forth.

I believe this is what LLR, PFGW, and prime95 all do for interim and final residues.

You can still do the no-mul-by-3's implementation. That is, N = 2^P-1, and you will calculate 3^(N+1) mod N. To print a standard iteration 1 residue, output 3^(top 3 bits of N+1)/3, etc. That is, you'll do divide-by-3 to output interim residues and divide-by-9 to get the final residue. This is what I'm doing in the next prime95 release -- testing it now.[/QUOTE]

When doing the mul-by-3 implementation, then for Mersenne numbers you should probably test 3^((N-1)/2) == -1, which is stronger than 3^(N-1) == 1.

R. Gerbicz 2017-09-21 07:56

[QUOTE=axn;468235]Because, your approach will not work for bases other than 2.[/QUOTE]
No and yes!
No, because my method works for any b, but yes because when we want to use this type
of ladder check we need to use b-th powermod calculation, this gives
that I will use (at least!) n*ceil(log2(b)) squaremod/mulmod, while
at a standard PRP test you are using (much) fewer, only n*log2(b) squaremods.

the core of the algorithm in b!=2:
[CODE]
u[i]=a^(b^i) mod N
d[t]=prod(i=0,t,u[L*i]) mod N
d[t]=d[t-1]*u[L*t] mod N
d[t]=a*d[t-1]^(b^L) mod N (to do the error check)
[/CODE]
With this a^(k*b^n)==u[n]^k mod N. (we delayed the k-th powermod).

New idea: you can save some work here: use say e=1000 and set B=b^e, L=1, but use error
checking after each 500 iteration. With such a huge B we reduce the overhead of the
b-th powermod to B-th powermod with (optimal or close to that) additon chain ([url]https://en.wikipedia.org/wiki/Addition-chain_exponentiation[/url])
And you need to do floor(n/e) iterations, after that do an additional (n mod e)
b-th powermod to get a^(b^n) mod N, and final k-th powermod.

But I'm not sure that even with this it is worth to do this error checking in the b!=2 test.

Ofcourse when b=2 we are already in an optimal case.
Multiple the final residue by a^(c-1) mod N if you want the Fermat residue if your number is N=k*b^n+c

[QUOTE=axn;468238]But that will not work in case of other bases, nor will it work in the case where there is a factor.[/QUOTE]

No. When you know a nontrivial d>1 factor of N=k*b^n+c then to test the primality of N/d, just
calculate a^N mod N what you would do in the d=1 case.
If (N/d) is Fermat prp number using base=a, then
a^(N/d)==a mod (N/d) this gives that
a^N==a^d mod (N/d),
so in the test get a^N mod N, then reduce it mod (N/d).
if it is not a^d mod (N/d), then N/d is not a Fermat prp number, the other direction is not true, because
a^N==a^d mod (N/d) doesn't give that
a^(N/d)==a mod (N/d), but it isn't a much weaker test.
Even in the thousand digits of range you wouldn't see a false hit in your life.
If this type of test passed then you can still do a true Fermat test for (N/d).

The above is a known trick from the standard Suyama test, what we use for example for
Fermat numbers, in the test, if we know the residue we get a fast prp test for N/d but that is also
weaker than a (single) Fermat test.

[QUOTE=preda;468239]When doing the mul-by-3 implementation, then for Mersenne numbers you should probably test 3^((N-1)/2) == -1, which is stronger than 3^(N-1) == 1.[/QUOTE]

For Mersenne prime(!) these are equivalent tests, for others I don't know even a single false hit, so at least in practice these are the same.

axn 2017-09-21 15:52

[QUOTE=R. Gerbicz;468245]No and yes!.[/QUOTE]

No, I don't mean that we can't implement your test. I mean that we can't do what George wants to do while also implementing this test. What George wants to do is, IIRC, to output residues (final & intermmediate) that corresponds to the computation 3^(N-1) for arbitrary N=(k*b^n+c)/f. Your test requires that we evaluate (3^k)^(b^n), and then equate it to 3^(f-c). These two sequences produce residues that are incompatible -- except for the case of f=1, b=2, where we can use George's trick to transform between the two.

R. Gerbicz 2017-09-21 16:28

[QUOTE=axn;468272]No, I don't mean that we can't implement your test. I mean that we can't do what George wants to do while also implementing this test. What George wants to do is, IIRC, to output residues (final & intermmediate) that corresponds to the computation 3^(N-1) for arbitrary N=(k*b^n+c)/f. Your test requires that we evaluate (3^k)^(b^n), and then equate it to 3^(f-c). These two sequences produce residues that are incompatible -- except for the case of f=1, b=2, where we can use George's trick to transform between the two.[/QUOTE]

Huge loss, but do you know that for f!=1 we have no (very) fast true primality test, like for example we still don't have a fast test for (10^n-1)/9 .

Btw my test could be even faster than George's test, if he really computes the Fermat residue for f>1, why? Because for b=2 I'm not doing mult by 2. *

[QUOTE=axn;468272]
Your test requires that we evaluate (3^k)^(b^n), and then equate it to 3^(f-c). [/QUOTE]
I'm prefering to write "a" instead of 3, because for Proth numbers we have a true primality test, and in those cases it is possible that 3 is a quadratic residue, so not usable.

* ps. the slowdown is much smaller if you do only a 64 bits mult once in 64 iterations. Maybe there could be other methods to reduce more this cost.

Prime95 2017-09-22 19:24

[QUOTE=R. Gerbicz;468245]
the core of the algorithm in b!=2:
[CODE]
u[i]=a^(b^i) mod N
[/CODE][/QUOTE]

Yes, I overlooked that this step requires the costly full mulmod rather than the easy mul-by-small-const. So, I concur that for b!=2 prime95 will not support Gerbicz error checks.

Prime95 2017-09-22 19:38

Here is my PRP residue proposal:

There are (at least) 5 PRP residue types:
1; Fermat PRP. N is PRP if a^(N-1) = 1 mod N
2: SPRP variant, N is PRP if a^((N-1)/2) = +/-1 mod N
3: Fermat PRP variant. N is PRP if a^(N+1) = a^2 mod N
4: SPRP variant. N is PRP if a^((N+1)/2) = +/-a mod N
5: Fermat PRP cofactor variant, N/d is PRP if a^(N-1) = a^d mod N/d

I encourage programs to return type 1 residues as that has been the standard for prime95, PFGW, LLR for many years.

I think Prime95 v29.4 will have options to support all of the above residue types for doing double-checks. PrimeNet will store residue type along with the residue -- need another JSON argument. When issuing PRP double-check assignments, the server will send the PRP base and desired residue type to the client. For first-time PRP assignments, the server will not specify a PRP base or residue type.

Interim residues are application dependent. However, it is strongly encouraged to report Fermat PRP residues using a standard binary powering algorithm. That is, iteration 0 is a, iteration 1 is a^(top 2 bits of N-1), etc.

The server will soon have about 50 type-4 residues from gpuOwl, which prime95 should be able to verify someday. While I'd prefer gpuOwl converts to type-1 residues, it is not a requirement.

The server also has several thousand type-1 PRP cofactor residues. I'm thinking of making type-5 the default for v29.4 as it has some advantages. 1) It allows Gerbicz error checking, 2) The residue remains the same as more and more factors are found for a Mersenne number - no double-check needed. 3) If a user decides to archive his PRP save files he can retest with a new known factor quite quickly.


All times are UTC. The time now is 21:16.

Powered by vBulletin® Version 3.8.11
Copyright ©2000 - 2021, Jelsoft Enterprises Ltd.