Williams' sequence 4*5^n-1 (A046865)

[Citrix]

Thank you for the excellent work. I am getting around 7000 kp/sec for the 3^16 sequence. That is 3.5 times faster than newpgen or jjsieve.

Is there a limit on the size of the k, your program can handle? I would like to test higher values like 3^32 and 3^64 etc.
Does your program have a prover code to submit primes found using your program.

[geoff]

Originally Posted by Citrix

Is there a limit on the size of the k, your program can handle? I would like to test higher values like 3^32 and 3^64 etc.

The current limit for k is 2^32-1. It is not a simple matter to increase that without affecting the sieve speed for general sequences.

It would be possible to write a seperate sieve routine that accepts k=a^b where a and b can each be up to 2^32-1. However there are other parts of the program that will limit the speed as b becomes large. In particular the prime sieve generates all primes and then checks to see if each one is of the correct form, this will become very ineffcient when b gets large, and in any case you will soon need to check factors larger than 2^62 which is the modular multiplication limit for the current code.

Have you tried one of these sieves?: GFNsieve, AthGFNsieve.

Does your program have a prover code to submit primes found using your program.

The entry is Srsieve.

[Citrix]

Originally Posted by geoff

The current limit for k is 2^32-1. It is not a simple matter to increase that without affecting the sieve speed for general sequences.

It would be possible to write a seperate sieve routine that accepts k=a^b where a and b can each be up to 2^32-1. However there are other parts of the program that will limit the speed as b becomes large. In particular the prime sieve generates all primes and then checks to see if each one is of the correct form, this will become very ineffcient when b gets large, and in any case you will soon need to check factors larger than 2^62 which is the modular multiplication limit for the current code.

Have you tried one of these sieves?: GFNsieve, AthGFNsieve.


The entry is Srsieve.

OK I will stick with this k and not go beyond 2^32

GFNsieve and AthGFNsieve only work with bases<2^32. For me the base is as large as 2^125,000. So I don't think they will work.

Thanks once again for modifying the sieve.

[Citrix]

I noticed that your program does not use the SPH algorithm. Do you have any plans of implementing it. It is quite straight forward. It should make your program several times faster.

It is described here
http://www.mersenneforum.org/showthread.php?t=4310
(See post by akruppa)

[geoff]

Originally Posted by Citrix

I noticed that your program does not use the SPH algorithm. Do you have any plans of implementing it. It is quite straight forward. It should make your program several times faster.

It is described here
http://www.mersenneforum.org/showthread.php?t=4310
(See post by akruppa)

Thanks for the links, I will have a look at this sometime. I have a feeling the mathematics behind it is more difficult than that behind baby-steps-giant-steps though, it might take me a while.

In srsieve version 0.3.14 I have made the code recognise numbers of the form A^(2^y)+B^(2^y) as generalised Fermats, and there is a bug fix to properly handle cases where the base is a square.

[geoff]

For the sequence 4*5^n-1, can anyone explain why all the factors end in 1 or 9 when n is odd?

I haven't tested many examples, but I think all the odd factors of k*5^n-1 end in 1 or 9 when n is odd and k is a square. Can anyone explain this?

[Citrix]

Originally Posted by geoff

Thanks for the links, I will have a look at this sometime. I have a feeling the mathematics behind it is more difficult than that behind baby-steps-giant-steps though, it might take me a while.

In srsieve version 0.3.14 I have made the code recognise numbers of the form A^(2^y)+B^(2^y) as generalised Fermats, and there is a bug fix to properly handle cases where the base is a square.

The alogrithm is the same as baby-steps-giant-steps algorithm only done in groups of factors of p-1, than based on sqrt of the range of n's.

so if p-1=2*a*b


then running time is sqrt(a)+sqrt(2)+sqrt(b) than sqrt(nmax-nmin)

Let me know if you have any questions.

[geoff]

Originally Posted by geoff

For the sequence 4*5^n-1, can anyone explain why all the factors end in 1 or 9 when n is odd?

I think I have found part of the answer to this: when n is odd, say n=2m+1, then 4*5^n-1 can be written as 5*x^2-y^2 where x=2*5^m and y=1, and according to a table I found in Guy, the primes with 5 as a quadratic residue are all -1,1 (mod 10).

I think the same idea (and the same table) can be used to speed up sieveing for the sequences 2*3^n-1, 5*6^n-1, 6*7^n-1, but I don't know how to work out the general case yet: Given k, what are the congruence classes of the primes with k as a quadratic residue?

[Citrix]

Originally Posted by geoff

The current limit for k is 2^32-1. It is not a simple matter to increase that without affecting the sieve speed for general sequences.

I think the limit is 2^31-1. I tried sieving k=2249392549 using your program it completely screws the k up, printing incorrect factors etc. If the limit is really 2^32 then I think this is a bug.

[geoff]

Originally Posted by Citrix

I think the limit is 2^31-1. I tried sieving k=2249392549 using your program it completely screws the k up, printing incorrect factors etc. If the limit is really 2^32 then I think this is a bug.

Thanks, this is a bug, I'll let you know when it is fixed.

[konrad127123]

Originally Posted by geoff

I think I have found part of the answer to this: when n is odd, say n=2m+1, then 4*5^n-1 can be written as 5*x^2-y^2 where x=2*5^m and y=1, and according to a table I found in Guy, the primes with 5 as a quadratic residue are all -1,1 (mod 10).

Correct :). I will try to explain where this table comes from.

Suppose p and q are odd primes. All odd primes are either 1 or 3 (mod 4). The law of Quadratic reciprocity says:

(A)If at least one of p and q is 1 (mod 4) then x^2=p (mod q) has a solution if and only if x^2=q (mod p) has a solution (so either both or neither of these congruences has a solution).

(B)If both p and q are 3 (mod 4) then exactly one of x^2=p (mod q) and x^2=q (mod p) has a solution (this case isn't relevant in the above)

Now let's set q=5. 5=1 (mod 4), so case A applies. Let p be the prime we're trial dividing by. So by case (A) x^2=5 (mod p) will only have a solution if and only if x^2=p (mod 5) has a solution. Substituting in x=0,1,2,3,4 in turn, we find that x^2 (mod 5) must be 0, 1 or 4. So if x^2=5 (mod p) has a solution then x^2 = p (mod 5) has a solution, so p=0,1 or 4 (mod 5).

So x^2=5 (mod p) can only have a solution if p=0,1 or 4 (mod 5). The only prime of the form p=0 (mod 5) is p=5 (which obviously won't divide any number of the form 4*5^n-1). So the only primes left are ones of the form p=1 or 4 (mod 5).

Now let's look at p (mod 10). If p is 1 or 4 (mod 5), then p is 1, 4, 6 or 9 (mod p). If p=4 or 6 (mod 10) then p is even and not equal to 2, so it is not prime. Thus p=1 or 9 (mod 10).

So if x^2 =5 (mod p) has a solution, then p=1 or 4 (mod p), so p=1 or 9 (mod 10), so we only need to consider such primes.

Originally Posted by geoff

I think the same idea (and the same table) can be used to speed up sieveing for the sequences 2*3^n-1, 5*6^n-1, 6*7^n-1, but I don't know how to work out the general case yet: Given k, what are the congruence classes of the primes with k as a quadratic residue?

If I remember correctly, a similar idea can be used for *all* sequences of the form k*b^n+/-c.Have a look at http://www.free-dc.org/forum/showthread.php?t=10262 (in particular look at post #6). It is written for base 2 Riesel's, but should be fairly easy to adapt to base 5. You would also need to read up on the Legendre symbol. (I can't remember off the top of my head if the Jacobi symbol is relevant here or not.)

If you have time, I recommend reading any good book on "Elementry Number Theory", which will cover such things.

Hope this helps :)