3-PRP faster than LL for GIMPS?

[bearnol]

import java.math.BigInteger;



public class Mersennef5 {

static BigInteger zero = new BigInteger("0");

static BigInteger one = new BigInteger("1");

static BigInteger two = new BigInteger("2");

static BigInteger three = new BigInteger("3");

static BigInteger nine = new BigInteger("9");

static BigInteger tenthousand = new BigInteger("10000");



public static void main (String args[])

throws java.io.IOException {



BigInteger n = new BigInteger("1");



BigInteger m = new BigInteger("0");



while (n.compareTo(tenthousand) < 0) {

if (isprime(n)) {

System.out.print(n.toString() + '\r');

m = mersenne(n);

if (fermat2(m))

System.out.println(n.toString());

}

n = n.add(one);

}

}



static BigInteger mersenne(BigInteger n) {

BigInteger mersenne = new BigInteger("1");

BigInteger i = new BigInteger("0");



while (i.compareTo(n) < 0) {

mersenne = mersenne.multiply(two);

i = i.add(one);

}



mersenne = mersenne.subtract(one);



return mersenne;

}

static boolean isprime(BigInteger n) {

return n.isProbablePrime(100);

}







static boolean fermat2(BigInteger n) {

boolean prime = true;





// workaround - apparent java bug in modPow - JGW

if (n.compareTo(one) == 0)

return false;



if (n.compareTo(two) == 0)

return true;



if (n.remainder(two).compareTo(zero) == 0)

return false;

// end workaround







if (power(three, (n.add(one)), n).compareTo(nine) == 0)

prime = true;

else

prime = false;



return prime;

}




static BigInteger power(BigInteger x, BigInteger n, BigInteger mod) {
// Knuth 4.6.3 - computes x^n modulo mod
BigInteger N = new BigInteger("0");
BigInteger Y = new BigInteger("0");
BigInteger Z = new BigInteger("0");

N=n;
Y=one;
Z=x;

while (true) {
if (N.remainder(two).compareTo(zero) > 0) {
N=N.divide(two);
Y=Z.multiply(Y).remainder(mod);
if (N.compareTo(zero) == 0)
return Y;
}
else {
N=N.divide(two);
}
Z=Z.multiply(Z).remainder(mod);
}
}

}

[R.D. Silverman]

Originally Posted by bearnol

import java.math.BigInteger;

<snip>

(1)

May I suggest that this forum is for the discussion of mathematics.
If you wish to discuss an algorithm, please give a mathematical description.
It is unreasonable to expect others to translate your code so that we may
understand what it does.

(2) Please explain why you think determining if a Mersenne prime candidate
is a 3-prp would be faster than a L-L test. Tell us what aspects of the
computation would be faster.

I will give a hint: It isn't. In fact, it would be slightly slower.
I will give a second hint: It won't prove primaility as L-L does.

[bearnol]

import java.math.BigInteger;



public class Mersenneu3 {

static BigInteger zero = new BigInteger("0");

static BigInteger one = new BigInteger("1");

static BigInteger two = new BigInteger("2");

static BigInteger tenthousand = new BigInteger("10000");



public static void main (String args[])

throws java.io.IOException {





BigInteger n = new BigInteger("1");



BigInteger m = new BigInteger("0");



while (n.compareTo(tenthousand) < 0) {

if (isprime(n)) {

System.out.print(n.toString() + '\r');

m = mersenne(n);

if (lucaslehmer(n, m))

System.out.println(n.toString());

}

n = n.add(one);

}

}



static BigInteger mersenne(BigInteger n) {

BigInteger mersenne = new BigInteger("1");

BigInteger i = new BigInteger("0");



while (i.compareTo(n) < 0) {

mersenne = mersenne.multiply(two);

i = i.add(one);

}



mersenne = mersenne.subtract(one);



return mersenne;

}



static boolean isprime(BigInteger n) {

return n.isProbablePrime(100);

}





static boolean lucaslehmer(BigInteger n, BigInteger m) {

BigInteger i = new BigInteger("1");

BigInteger S = new BigInteger("4");



for (i = one; i.compareTo(n.subtract(one)) < 0; i = i.add(one))

S = ((S.multiply(S)).subtract(two)).remainder(m);



if (S.compareTo(zero) == 0)

return true;

else

return false;

}







}

[bearnol]

I posted the two separate competing sections of code first, on its own, and in its entirety, in case anyone fancied repeating my experiment unchanged. Clearly there are other issues (eg FFT multiplication), but maybe these apply equally.
My math reasoning goes like this:
First, you will notice I use the "Russian Peasant" method for evaluating modPow in the 3-PRP code. Secondly I actually evaluate 3^(n+1) [mod n] (rather than 3^n) and check for equality with nine (rather than 3). This hopefully cunning alteration means the modPow, in the case of a Mersenne number, is being raised to an exact power of two, which should further speed it up.
Obviously a 3-PRP is not as definitive in _proving_ primality, but _is_ as definitive in _disproving_ it, which is what we want for GIMPS 99% of the time.
J

[R.D. Silverman]

Originally Posted by bearnol

I posted the two separate competing sections of code first, on its own, and in its entirety, in case anyone fancied repeating my experiment unchanged. Clearly there are other issues (eg FFT multiplication), but maybe these apply equally.
My math reasoning goes like this:
First, you will notice I use the "Russian Peasant" method for evaluating modPow in the 3-PRP code. Secondly I actually evaluate 3^(n+1) [mod n] (rather than 3^n) and check for equality with nine (rather than 3). This hopefully cunning alteration means the modPow, in the case of a Mersenne number, is being raised to an exact power of two, which should further speed it up.
Obviously a 3-PRP is not as definitive in _proving_ primality, but _is_ as definitive in _disproving_ it, which is what we want for GIMPS 99% of the time.
J

Just what we need. Another crank.

I requested that the poster describe the algorithm, rather than posting code.
So what does he do? He posts yet more code. And people wonder why I lose
patience.

I suggest to the poster that he show EXACTLY how his powering algorithm
reduces the number of multiplications. [Hint: it doesn't; I told him this
in private email, but of course this was ignored and followed by another
public post]

His trick of computing 3^(n+1), rather than 3^n is indeed useful since
every iteration is now just a squaring. But it STILL REQUIRES n
multiplications! [squarings]

Go read the section on "addition chains" in Knuth Vol 2. One can not
do better than n squarings. log_2[N] is a lower bound on the length of
any addition chain. When N = 2^n, we get a chain of length n.

[bearnol]

1) With the 3-PRP, we no longer need the subtraction of 2 at each step (as per LL), so maybe we gain here?
2) I've done the actual test [at least in Java, up to n=10000], and it seems to work, so I was looking for explanation, as much as anything else...
J

[R.D. Silverman]

Originally Posted by bearnol

1) With the 3-PRP, we no longer need the subtraction of 2 at each step (as per LL), so maybe we gain here?
2) I've done the actual test [at least in Java, up to n=10000], and it seems to work, so I was looking for explanation, as much as anything else...
J

(1) The time to do the subtraction is lost in the noise. OTOH, the "russian
peasant" powering algorithm has a fair amount of non-constant overhead per
iteration beyond that of the single precision subtraction in L-L.

(2) "seems to work". What seems to work? Of course one can do a base 3
prp test. And If you wanted an explanation for something, then why did you post CODE???

If you are asking whether 3 is ever a FALSE witness to 2^n-1, the answer
is certainly yes, but the occurrence will be VERY VERY rare. On rough
probabilistic grounds I would expect #{x < N | x = 2^n-1, x composite, 3 is a false witness} to be no more than loglog N, i.e. false witnesses should
be at least as rare as Mersenne primes themselves and certainly as hard to
find.

[ppo]

Bearnol, can you please put comments in your code, if you want us to understand it ? (no, JAVA is not a self documenting language).
beside that, the entire thread should probably go to the software forum, since it is dealing mostly with implementation details.
thanks. Pierpaolo

[flava]

IMHO, a 3-PRP test (or any k-PRP) and the LL test for a Mersenne M(n) have the same complexity O(n) if they use the same method for multiplication (squaring). The lack of substraction for each step has almost zero influence on the total time, a substraction takes several orders of magnitude less cycles than a multiplication. Make the following test:
- supose you want to test M(n)
- calculate M(n) - 2 (n-1 times) and measure the total time
=> this is the maximum gain you can expect

[bearnol]

...so in fact (I think we're now all agreed :-) ) PRP would prove (slightly) faster than LL [as it would not need the subtractions, but otherwise be similar]. This is what I have observed in practice. Who's going to tell George Woltman?!

[akruppa]

No one will have to tell him, he reads the forum himself. The PRP idea has been suggested many times before, i.e. on the Mersenne mailing list before the forum even existed.

The 3-PRP is NOT faster than a LL. The time for the subtraction at each step is completely negligible. And we'll need the LL test anyway so we can actually prove primality of Mersenne numbers. Why should we add another algorithm that is exactly as fast as LL but inferior in that it is not sufficient for a primality proof?

Alex