## Re: Generating (small?) 8-16K digit prime numbers

Expand Messages
• Sirs, In reply, the mpz_probab_prime_p() function performs some trial division, a simple fermat test, then a miller-rabin test, whilst the mpz_frobenius test
Message 1 of 9 , Jun 16 3:14 PM
Sirs,

In reply, the mpz_probab_prime_p() function performs some trial
division, a simple fermat test, then a miller-rabin test, whilst the
mpz_frobenius test performs a fast test for a perfect power then on to
it's main test, so I believe I have covered the basics.

As for the wheel method, I can only assume that trial division up to
the number if binary digits / 30 of the said number (which I
understand is what GMP does) is sufficient to rule out a lot of
possible primes.

I humbly agree to the rather bad choice of the last digit (possibly
including 5 -- oops! although this will be ruled out very early with
trial division).

Another poster suggested I try using something like proth, however a)
I am looking for primes with no specialized form (2^k+/-n etc) and b)
the intention is to have a simple command line program which I can run
in the background in *nix/windows.

Thanks anyway Alan (and others) for your assistance.

--- In primenumbers@yahoogroups.com, Alan Eliasen <eliasen@m...> wrote:
>
> Alan McFarlane wrote:
> > Sirs
> >
> > I've been looking around various sources for a method of quickly
> > generating primes in the 8,000 to 16,000 odd (decimal) digit range
and
> > am finding some difficulty.
>
> Hi Alan (good name!)
>
> These tests become quite expensive at the 8,000 digit level. Can I
> offer a few optimizations?
>
> * Try some trial divisions first with small primes. (You can skip 2
> given your method of choosing last number.) These will be far less
> expensive than doing a full probable prime test or strong probable prime
> test on these numbers. You'll be able to very quickly reject a large
> percentage of numbers.
>
> * Your method of choosing the last digit could be slightly improved.
> You can eliminate "5" as a possible last digit. Maybe simply choosing
> from the final digits [1379] would help. But trial divisions by small
> primes would make this largely unnecessary.
>
> * Look into the "wheel" method of eliminating composite numbers so
> you don't generate them as possibilities in the first place. This one
> is rather simple. Here's how to understand what's happening. Take out
> a sheet of graph paper and divide it into 30 (hint: 30 is 2*3*5)
> columns. It's like a Sieve of Eratosthenes. The first row represents
> the numbers 0-29, the second row represents the numbers 30-59, etc. Put
> an "x" through all the squares that are divisible by smaller primes.
> That is, 2,4,6,... then 3,6,9..., then 5,10,15...
>
> You'll see a pattern in the squares that are left. The columns that
> are always crossed out can't generate primes. So, you can generate
> numbers that are more probably prime by:
>
> rand[x] * 30 + (one of) [1,7,11,13,17,19,23,29]
>
> This is largely equivalent to trial-factoring by 2,3, and 5, but can
> potentially save even more time. Hope this helps.
>
> --
> Alan Eliasen | "It is not enough to do your best;
> eliasen@m... | you must know what to do and THEN
> http://futureboy.homeip.net/ | do your best." -- W. Edwards Deming
>
>
> > I've been using a miller-rabin test and a frobenius test (both with 5
> > repetitions) however I have discovered that my program is *extremely*
> > slow. It took several days just to generate one 8,000 digit prime!
> >
> > [code snippet using MinGW with GMP - error checking removed]
> >
> > void generate_prime( mpz_t rop, unsigned digits )
> > {
> > static char charset[] = "0123456789";
> > char* buffer = (char*)malloc(digits + 1);
> > unsigned d;
> >
> > do
> > {
> > do
> > {
> > /* first digit - 1..9 inclusive */
> > do
> > {
> > buffer[0] = charset[(rand() % 9) + 1];
> > }
> > while (buffer[0] == '0');
> >
> > /* middle digits - 0..9 inclusive */
> > for (d = 1; d < digits - 1; d++)
> > buffer[d] = charset[rand() % 10];
> >
> > /* last digit */
> > buffer[digits - 1] = charset[((rand() % 5) * 2) + 1];
> >
> > /* trailing '\0' */
> > buffer[digits] = '\0';
> >
> > /* set 'rop' */
> > mpz_set_str(rop, buffer, 10);
> > }
> > while (mpz_probab_prime_p(rop, 5) == 0);
> > }
> > while (mpz_frobenius(rop, 5) == 0);
> >
> > free(buffer);
> > }
> >
> > [/code]
> >
> > This code, whilst it works in-as-much, as I have verified various
sizes
> > of primes up to 2000 digits successfully, is extremely slow.
> >
> > I can only assume that I am being very naive in choosing such an
> > algorithm. Perhaps someone could suggest a) a different algorithm, b)
> > some appropriate reading material (assuming college level mathematics
> > albeit a long time ago!) or c) some links to where I may learn more
> >
> >
> > --
> > Alan McFarlane
• ... Yes, it does. Primes are sparse up at 8000 digits (if you just purely chose random numbers, only about 1 in 18000 would be prime.) I think it ll just
Message 2 of 9 , Jun 16 4:27 PM
mcfarlane_alan wrote:
> Sirs,
>
> In reply, the mpz_probab_prime_p() function performs some trial
> division, a simple fermat test, then a miller-rabin test, whilst the
> mpz_frobenius test performs a fast test for a perfect power then on to
> it's main test, so I believe I have covered the basics.

Yes, it does. Primes are sparse up at 8000 digits (if you just
purely chose random numbers, only about 1 in 18000 would be prime.) I
think it'll just take a long time, assuming that GMP's implementations
are sensible, but there are a few things you can test:

Have you done a test to see how fast the Rabin-Miller test is
compared to the Frobenius test? Rabin-Miller is generally one of the
fastest algorithms. I don't know exactly how the Frobenius test scales
in comparison, both practically and asymptotically.

You could maybe do strong pseudoprime tests base 2 and 3, and then a
Lucas test. This combination is often known as the Baillie-PSW which
has no known failures (but it's possible that there are some.)

You're doing only 5 rounds of Rabin-Miller which seems like a
reasonable number given your bit sizes. While the worst case of
Rabin-Miller is usually cited as giving a probability of 4^-n of
erroneously declaring a composite number probably prime, when n are the
number of rounds, but the probability of failure usually decreases
greatly at large bit sizes. Paper 3 below gives estimates for
probabilities, but I don't recall that they went anywhere as high as
8000 digits. The other papers are interesting resources for discussing
cases where Rabin-Miller may fail.

1. W. R. Alford, A. Granville and C. Pomerance, On the difficulty of
finding reliable witnesses, Algorithmic Number Theory, pp. 1-16, Lecture
Notes in Computer Science, vol. 877, Springer-Verlag, Berlin, 1994. MR
96d:11136

2. François Arnault, Constructing Carmichael numbers which are strong
pseudoprimes to several bases, Journal of Symbolic Computation, v.20
n.2, p.151-161, Aug. 1995

3. I. Damgård, P. Landrock, and C. Pomerance, Average case estimates
for the strong probable prime test, Math. Comp. 61 (1993), 177-194. MR
94b:11124

4. G. Jaeschke, On strong pseudoprimes to several bases, Math. Comp. 61
(1993), 915-926. MR 94d:11004

5. Louis Monier, Evaluation and comparison of two efficient
probabilistic primality testing algorithms , Theoretical Computer
Science 12 (1980), 97-108. MR 82a:68078

In any case, let us know what you find. This is a common question
for lots of implementers.

--
Alan Eliasen | "It is not enough to do your best;
eliasen@... | you must know what to do and THEN
http://futureboy.homeip.net/ | do your best." -- W. Edwards Deming

> As for the wheel method, I can only assume that trial division up to
> the number if binary digits / 30 of the said number (which I
> understand is what GMP does) is sufficient to rule out a lot of
> possible primes.
>
> I humbly agree to the rather bad choice of the last digit (possibly
> including 5 -- oops! although this will be ruled out very early with
> trial division).
>
> Another poster suggested I try using something like proth, however a)
> I am looking for primes with no specialized form (2^k+/-n etc) and b)
> the intention is to have a simple command line program which I can run
> in the background in *nix/windows.
>
> Thanks anyway Alan (and others) for your assistance.
>
> --- In primenumbers@yahoogroups.com, Alan Eliasen <eliasen@m...> wrote:
>
>>Alan McFarlane wrote:
>>
>>>Sirs
>>>
>>>I've been looking around various sources for a method of quickly
>>>generating primes in the 8,000 to 16,000 odd (decimal) digit range
>
> and
>
>>>am finding some difficulty.
>>
>> Hi Alan (good name!)
>>
>> These tests become quite expensive at the 8,000 digit level. Can I
>>offer a few optimizations?
>>
>> * Try some trial divisions first with small primes. (You can skip 2
>>given your method of choosing last number.) These will be far less
>>expensive than doing a full probable prime test or strong probable prime
>>test on these numbers. You'll be able to very quickly reject a large
>>percentage of numbers.
>>
>> * Your method of choosing the last digit could be slightly improved.
>> You can eliminate "5" as a possible last digit. Maybe simply choosing
>>from the final digits [1379] would help. But trial divisions by small
>>primes would make this largely unnecessary.
>>
>> * Look into the "wheel" method of eliminating composite numbers so
>>you don't generate them as possibilities in the first place. This one
>>is rather simple. Here's how to understand what's happening. Take out
>>a sheet of graph paper and divide it into 30 (hint: 30 is 2*3*5)
>>columns. It's like a Sieve of Eratosthenes. The first row represents
>>the numbers 0-29, the second row represents the numbers 30-59, etc. Put
>>an "x" through all the squares that are divisible by smaller primes.
>>That is, 2,4,6,... then 3,6,9..., then 5,10,15...
>>
>> You'll see a pattern in the squares that are left. The columns that
>>are always crossed out can't generate primes. So, you can generate
>>numbers that are more probably prime by:
>>
>> rand[x] * 30 + (one of) [1,7,11,13,17,19,23,29]
>>
>> This is largely equivalent to trial-factoring by 2,3, and 5, but can
>>potentially save even more time. Hope this helps.
>>
>>--
>> Alan Eliasen | "It is not enough to do your best;
>> eliasen@m... | you must know what to do and THEN
>> http://futureboy.homeip.net/ | do your best." -- W. Edwards Deming
>>
>>
>>
>>>I've been using a miller-rabin test and a frobenius test (both with 5
>>>repetitions) however I have discovered that my program is *extremely*
>>>slow. It took several days just to generate one 8,000 digit prime!
>>>
>>>[code snippet using MinGW with GMP - error checking removed]
>>>
>>>void generate_prime( mpz_t rop, unsigned digits )
>>>{
>>> static char charset[] = "0123456789";
>>> char* buffer = (char*)malloc(digits + 1);
>>> unsigned d;
>>>
>>> do
>>> {
>>> do
>>> {
>>> /* first digit - 1..9 inclusive */
>>> do
>>> {
>>> buffer[0] = charset[(rand() % 9) + 1];
>>> }
>>> while (buffer[0] == '0');
>>>
>>> /* middle digits - 0..9 inclusive */
>>> for (d = 1; d < digits - 1; d++)
>>> buffer[d] = charset[rand() % 10];
>>>
>>> /* last digit */
>>> buffer[digits - 1] = charset[((rand() % 5) * 2) + 1];
>>>
>>> /* trailing '\0' */
>>> buffer[digits] = '\0';
>>>
>>> /* set 'rop' */
>>> mpz_set_str(rop, buffer, 10);
>>> }
>>> while (mpz_probab_prime_p(rop, 5) == 0);
>>> }
>>> while (mpz_frobenius(rop, 5) == 0);
>>>
>>> free(buffer);
>>>}
>>>
>>>[/code]
>>>
>>>This code, whilst it works in-as-much, as I have verified various
>
> sizes
>
>>>of primes up to 2000 digits successfully, is extremely slow.
>>>
>>>I can only assume that I am being very naive in choosing such an
>>>algorithm. Perhaps someone could suggest a) a different algorithm, b)
>>>some appropriate reading material (assuming college level mathematics
>>>
>>>
>>>--
>>>Alan McFarlane
• ... Both your criteria are satisfied by PrimeFormGW. This is a very efficient and sophisticated program, which also takes care of preliminary trial division in
Message 3 of 9 , Jun 17 2:26 AM
In an email dated 16/6/2005 11:14:36 pm GMT Daylight time, "mcfarlane_alan" <alan.mcfarlane@...> writes:

>Another poster suggested I try using something like proth, however a)
>I am looking for primes with no specialized form (2^k+/-n etc) and b)
>the intention is to have a simple command line program which I can run
>in the background in *nix/windows.
>

Both your criteria are satisfied by PrimeFormGW.

This is a very efficient and sophisticated program, which also takes care of preliminary trial division in an intelligent and configurable way.

It will never miss a prime, and (for a given choice of base) will only throw up a few pseudoprimes. These can then be further investigated by any one of several methods, including your favoured Miller-Rabin.

-Mike Oakes
• From: Alan McFarlane ... Firstly, if possible, sieve your candidates. Then use PFGW, A.K.A. OpenPFGW, ne e primeform. The primeform
Message 4 of 9 , Jun 17 3:18 AM
From: Alan McFarlane <alan.mcfarlane@...>
> Sirs
>
> I've been looking around various sources for a method of quickly
> generating primes in the 8,000 to 16,000 odd (decimal) digit range and
> am finding some difficulty.

Firstly, if possible, sieve your candidates.
Then use PFGW, A.K.A. OpenPFGW, ne\'e primeform.
The primeform yahoogroup is a sister group to this one, exes for windows and

> I've been using a miller-rabin test and a frobenius test (both with 5
> repetitions)

That's 4 too many repetitions:
http://primes.utm.edu/notes/prp_prob.html

5 repetitions doesn't noticably increase the probability of the number being
prime - 99.<800-nines>% vs. 99.<4000-nines>% (alethic issues notwithstanding).

Similarly, it doesn't _at all_ increase your knowledge whether the number is
actually prime, which starts at zero and ends up being zero.

Phil

() ASCII ribbon campaign () Hopeless ribbon campaign
/\ against HTML mail /\ against gratuitous bloodshed

[stolen with permission from Daniel B. Cristofani]

__________________________________________________
Do You Yahoo!?
Tired of spam? Yahoo! Mail has the best spam protection around
http://mail.yahoo.com
• Sirs, ... Ah, good point. Perhaps I should start with a random number, and incremently test from there onwards. This would mean that (in theory) I would hit a
Message 5 of 9 , Jun 19 12:46 AM
Sirs,

>In reply, the mpz_probab_prime_p() function performs some trial
>>division, a simple fermat test, then a miller-rabin test, whilst the
>>mpz_frobenius test performs a fast test for a perfect power then on to
>>it's main test, so I believe I have covered the basics.

>Yes, it does. Primes are sparse up at 8000 digits (if you just
>purely chose random numbers, only about 1 in 18000 would be prime.) I
>think it'll just take a long time, assuming that GMP's implementations
>are sensible, but there are a few things you can test:

Ah, good point. Perhaps I should start with a random number, and
incremently test from there onwards. This would mean that (in theory) I
would hit a prime in ~18000 tests, whereas testing purely random numbers
may result in more tests. (Not sure on this one though).

>Have you done a test to see how fast the Rabin-Miller test is
>compared to the Frobenius test? Rabin-Miller is generally one of the
>fastest algorithms. I don't know exactly how the Frobenius test scales
>in comparison, both practically and asymptotically.

I believe a frobenius test is roughly 3 time slower than a miller-rabin
test, although the frobenius test has a much better failure rate - cited
as being less than 1/1710.

>>Another poster suggested I try using something like proth, however a)
>>I am looking for primes with no specialized form (2^k+/-n etc) and b)
>>the intention is to have a simple command line program which I can run
>>in the background in *nix/windows.

>Both your criteria are satisfied by PrimeFormGW.

>This is a very efficient and sophisticated program, which also takes
>care of preliminary trial division in an intelligent and configurable
>way.

>It will never miss a prime, and (for a given choice of base) will only
>throw up a few pseudoprimes. These can then be further investigated by
>any one of several methods, including your favoured Miller-Rabin.

Thanks Mike, I'll have a good look at this.

>>I've been using a miller-rabin test and a frobenius test (both with 5
>>repetitions)

>That's 4 too many repetitions:
>http://primes.utm.edu/notes/prp_prob.html

Woah! Nice one Phil. I'm sure Knuth, Donald suggested that 25 repetions

>A few days ago, I've posted a commentary at "The GMP_Discussion List"
>about the starting method I'm intending to use which is equivalent to
>TRIAL DIVISIONS for all 203,280,220 primes <= 2^32-5 = 4,294,967,291
>(check it at
>http://swox.com/list-archives/gmp-discuss/2005-June/001692.html ). In
>the worst case, my actual routines require about one and half hour to
>test rare numbers that are multiple of only primes next to the above
>cited limit, however, for the majority (which has some small prime
>factors), it takes an average of only a few seconds !!!

Hmm, I can easily perform a trial division to 2^64 odd, however I tend
to run trial division to around 2^16 (for speed) then perform the
light-weight SPRP tests and then onto the heavy ones (such as frobenius).

Thanks again for you help people, I'll let you know what my results are
after a brief holiday.

--
Alan
• I would recommend looking at the GPG and PGP code for generating random numbers. I believe you can pick how large the numbers are in the library call. They use
Message 6 of 9 , Jun 19 5:38 AM
I would recommend looking at the GPG and PGP code for generating random numbers. I
believe you can pick how large the numbers are in the library call. They use Fermat for
the first phase which is very fast, and then something else in the second phase. It has
been a few years since I looked it up, but I was impressed when I ran it and examined the
code.

Ian Smith

__________________________________________________
Do You Yahoo!?
Tired of spam? Yahoo! Mail has the best spam protection around
http://mail.yahoo.com
Your message has been successfully submitted and would be delivered to recipients shortly.