Post Reply 
[HP35s] Program for prime number (Wheel Sieve and Miller-Rabin)
02-16-2019, 10:45 AM
Post: #21
RE: [HP35s] Program for prime number (brut force)
(02-13-2019 11:53 AM)Thomas Klemm Wrote:  You may want to have a look at (35S) Number Theory Library.
There are the following functions:

Quote:PRIME? 148
A -> 0 OR 1 as A is composite or probably prime.

RABIN? 571
A, B -> 0 or 1, 1 if A is a strong pseudo-prime base B, 0 if not.

Though I haven't tried but if PRIME? uses RABIN? (which I assume implements the Miller-Rabin Primality Test) it might be considerably faster for bigger numbers.

I've implemented this for the HP-48GX:
(48G) Miller-Rabin Primality Test

This article gives a bit more details:
Miller-Rabin Primality Test for the HP-48

Cheers
Thomas

Yes, PRIME? calls RABIN? to one random base after first checking for small factors of A.
Find all posts by this user
Quote this message in a reply
02-16-2019, 11:16 AM
Post: #22
RE: [HP35s] Program for prime number (brut force)
The PRIME? programme in

http://www.hpmuseum.org/forum/thread-4236.html

for input

999 999 999 989

returns

1

ie entry is prime, in 36 seconds.
Find all posts by this user
Quote this message in a reply
02-16-2019, 11:36 AM
Post: #23
RE: [HP35s] Program for prime number (brut force)
(02-11-2019 05:16 PM)fred_76 Wrote:  
PHP Code:
Digits        Prime Number           Total Time
12 digits      999 999 999 989       02
:32:18 

(02-16-2019 11:16 AM)Gerald H Wrote:  for input 999 999 999 989 returns 1 ie entry is prime, in 36 seconds.

Impressive: That's 9138 ÷ 36 = 253.83 times faster.
Find all posts by this user
Quote this message in a reply
02-16-2019, 11:48 AM
Post: #24
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 11:36 AM)Thomas Klemm Wrote:  
(02-11-2019 05:16 PM)fred_76 Wrote:  
PHP Code:
Digits        Prime Number           Total Time
12 digits      999 999 999 989       02
:32:18 

(02-16-2019 11:16 AM)Gerald H Wrote:  for input 999 999 999 989 returns 1 ie entry is prime, in 36 seconds.

Impressive: That's 9138 ÷ 36 = 253.83 times faster.

True ! But I’m pretty sure it can still be optimized. I can’t do that right now. See you in 2-3 weeks...
Find all posts by this user
Quote this message in a reply
02-16-2019, 11:57 AM
Post: #25
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 10:45 AM)Gerald H Wrote:  Yes, PRIME? calls RABIN? to one random base after first checking for small factors of A.

Thus in case of 1 we only know that it's a strong probable prime.
But we could be sure by checking the bases 2, 13, 23, and 1662803.

Cheers
Thomas
Find all posts by this user
Quote this message in a reply
02-16-2019, 12:54 PM
Post: #26
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 11:57 AM)Thomas Klemm Wrote:  
(02-16-2019 10:45 AM)Gerald H Wrote:  Yes, PRIME? calls RABIN? to one random base after first checking for small factors of A.

Thus in case of 1 we only know that it's a strong probable prime.
But we could be sure by checking the bases 2, 13, 23, and 1662803.

Cheers
Thomas

Please report a number falsely returned as prime.
Find all posts by this user
Quote this message in a reply
02-16-2019, 02:12 PM
Post: #27
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 12:54 PM)Gerald H Wrote:  Please report a number falsely returned as prime.

Example:
Quote:Suppose we wish to determine if n = 221 is prime.
We randomly select a number a such that 1 < a < n - 1, say a = 174.

Since 220 ≡ −1 mod n, either 221 is prime, or 174 is a strong liar for 221.

But 221 = 13×17, so it's not a prime.

Quote:It can be shown that for any odd composite n, at least 3/4 of the bases a are witnesses for the compositeness of n.

Thus chances are that a composite number is falsely considered a prime if we run the test only once.
I would expect this to happen if you run the program multiple times with the same composite number.

Cheers
Thomas
Find all posts by this user
Quote this message in a reply
02-16-2019, 02:22 PM
Post: #28
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 12:54 PM)Gerald H Wrote:  
(02-16-2019 11:57 AM)Thomas Klemm Wrote:  Thus in case of 1 we only know that it's a strong probable prime.
But we could be sure by checking the bases 2, 13, 23, and 1662803.

Cheers
Thomas

Please report a number falsely returned as prime.
It seems hard because the algorithm uses a random test. Therefore on number may be missed one time but not the other time. By testing 2/13/23/1662803, you are sure that (up to 12 digits) your number is prime or composite.
Find all posts by this user
Quote this message in a reply
02-16-2019, 02:29 PM (This post was last modified: 02-16-2019 03:02 PM by Albert Chan.)
Post: #29
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 10:45 AM)Gerald H Wrote:  Yes, PRIME? calls RABIN? to one random base after first checking for small factors of A.

Just checking, when you say random base, you meant 2 to A-2 ?

Technically, this is a composite test. Any number that failed is *guaranteed* composite.

The reverse is not always true: strong pseudoprimes

Thus, we use pre-selected set of bases, to cover each other's "holes", to prove primality.

Baillie-PSW primality test (SPRP + Lucas test) covered the "holes" even better.
There is an award of $30 to show an example of pseudoprime.
Find all posts by this user
Quote this message in a reply
02-16-2019, 05:39 PM (This post was last modified: 02-16-2019 05:42 PM by Gerald H.)
Post: #30
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 02:12 PM)Thomas Klemm Wrote:  
(02-16-2019 12:54 PM)Gerald H Wrote:  Please report a number falsely returned as prime.

Example:
Quote:Suppose we wish to determine if n = 221 is prime.
We randomly select a number a such that 1 < a < n - 1, say a = 174.

Since 220 ≡ −1 mod n, either 221 is prime, or 174 is a strong liar for 221.

But 221 = 13×17, so it's not a prime.

Quote:It can be shown that for any odd composite n, at least 3/4 of the bases a are witnesses for the compositeness of n.

Thus chances are that a composite number is falsely considered a prime if we run the test only once.
I would expect this to happen if you run the program multiple times with the same composite number.

Cheers
Thomas

"Chances are..." - No.

"At least 3/4 of bases are witnesses..." - Yes, but how many are in fact witnesses? The "at least" allows the likelihood that many more than 3/4 are witnesses, & that is in fact generally the case.

Set at random

N := 803189 * 485909

the product of two primes, how many bases are actually witnesses? Try running PRIME?, I will be very surprised if the 35s is fast enough for you to repeat the test until a 1 appears.

I use one test as it's reliable enough - Certainty is expensive, if you're really not sure try running the programme again.

For larger numbers, say 2^555 range, one test gives certainty.
Find all posts by this user
Quote this message in a reply
02-16-2019, 05:43 PM
Post: #31
RE: [HP35s] Program for prime number (brut force)
Additionally, PRIME? removes 221 via the small factors test.
Find all posts by this user
Quote this message in a reply
02-16-2019, 05:57 PM
Post: #32
RE: [HP35s] Program for prime number (brut force)
This is a list of composite numbers \(n < 10000\) with the amount of strong liars \(k\) and the ratio \(r=\frac{n-3}{k}\):

Code:
   r    n    k
   4.21 1891 448
   5.00 8911 1780
   5.57 2701 484
   6.18 6533 1056
   6.20 5461 880
   6.41 1541 240
   8.22 8321 1012
   8.33 4033 484
   8.52 2047 240
   8.65 1387 160

Thus I suggest to run PRIME? multiple time with 1891 and see how long it takes until you hit a strong liar.

Cheers
Thomas
Find all posts by this user
Quote this message in a reply
02-16-2019, 05:59 PM
Post: #33
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 05:57 PM)Thomas Klemm Wrote:  This is a list of composite numbers \(n < 10000\) with the amount of strong liars \(k\) and the ratio \(r=\frac{n-3}{k}\):

Code:
   r    n    k
   4.21 1891 448
   5.00 8911 1780
   5.57 2701 484
   6.18 6533 1056
   6.20 5461 880
   6.41 1541 240
   8.22 8321 1012
   8.33 4033 484
   8.52 2047 240
   8.65 1387 160

Thus I suggest to run PRIME? multiple time with 1891 and see how long it takes until you hit a strong liar.

Cheers
Thomas

1891 will not get that far, declared composite by the small factors test.
Find all posts by this user
Quote this message in a reply
02-16-2019, 06:07 PM
Post: #34
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 05:39 PM)Gerald H Wrote:  1891 will not get that far, declared composite by the small factors test.

So what's the biggest of your small factors that you test?
Find all posts by this user
Quote this message in a reply
02-16-2019, 06:18 PM
Post: #35
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 06:07 PM)Thomas Klemm Wrote:  
(02-16-2019 05:39 PM)Gerald H Wrote:  1891 will not get that far, declared composite by the small factors test.

So what's the biggest of your small factors that you test?

Have you actually had a look at or used the programme on the 35s?

Should you inspect the programme you'll find the largest small factor that is tested for is

999999
Find all posts by this user
Quote this message in a reply
02-16-2019, 07:09 PM
Post: #36
RE: [HP35s] Program for prime number (brut force)
"Set at random

N := 803189 * 485909

the product of two primes, how many bases are actually witnesses? Try running PRIME?, I will be very surprised if the 35s is fast enough for you to repeat the test until a 1 appears."

The use of "witness" above is correct, most bases return 0 for input N & are thus "witnesses".

Please post any comments in this thread for the benefit of all.
Find all posts by this user
Quote this message in a reply
02-16-2019, 07:57 PM (This post was last modified: 02-16-2019 08:52 PM by Albert Chan.)
Post: #37
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 07:09 PM)Gerald H Wrote:  "Set at random

N := 803189 * 485909

the product of two primes, how many bases are actually witnesses? Try running PRIME?, I will be very surprised if the 35s is fast enough for you to repeat the test until a 1 appears."

Re-reading your quote, I was mistaken.
I read a bit too fast, and think you are saying until 1 [witness] appears.

Changing the word witness to non-witness made the statement un-ambiguous.
"Until 1 pseduoprime appears" also work.


(02-16-2019 06:18 PM)Gerald H Wrote:  Should you inspect the programme you'll find the largest small factor that is tested for is

999999

If the code checked this far, why the need for SPRP test ?
Any 12 digits integer, with small factor upto 999999 checked, already proved primality.
Find all posts by this user
Quote this message in a reply
02-16-2019, 08:39 PM
Post: #38
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 05:39 PM)Gerald H Wrote:  Set at random

N := 803189 * 485909

the product of two primes, how many bases are actually witnesses?

It could very well be that this number has no strong liars.

(02-16-2019 06:18 PM)Gerald H Wrote:  Have you actually had a look at or used the programme on the 35s?

I had a look at the program but wasn't able to follow.
Even if the batteries of my HP-35s weren't dead I doubt I would enter an 850 line program just to figure that out by myself.

Quote:Should you inspect the programme you'll find the largest small factor that is tested for is

999999

Not sure if I understood that correctly but I checked products \(n = p \cdot q\) of primes \(p, q > 1000\).
So they shouldn't be detected by testing small factors.
Code:
      r      n       k
      417.11 1102837 2644
     1012.01 1018081 1006
     1016.01 1026169 1010
     1022.01 1038361 1016
     1024.01 1042441 1018
     1034.01 1062961 1028
     1036.01 1067089 1030
     1042.01 1079521 1036
     1052.01 1100401 1046
     1054.01 1104601 1048
     1064.01 1125721 1058
     1066.01 1129969 1060
     1072.01 1142761 1066
     1090.01 1181569 1084
     1094.01 1190281 1088
     1096.01 1194649 1090
     1100.01 1203409 1094
     1205.06 1060459 880
     1305.39 1148743 880
     2395.24 1073071 448
     5317.83 1042297 196
     6703.52 1072567 160
     7319.45 1083281 148

Chances might be small with \(1102837=1009×1093\) to hit a strong liar but they are not 0.
I was only aware of the upper bound of the probability \(p<\frac{1}{4}\) but not of the exact distribution.

Cheers
Thomas
Find all posts by this user
Quote this message in a reply
02-16-2019, 10:32 PM (This post was last modified: 02-16-2019 10:40 PM by Albert Chan.)
Post: #39
RE: [HP35s] Program for prime number (brut force)
(02-16-2019 05:57 PM)Thomas Klemm Wrote:  
Code:
   r    n    k
   4.21 1891 448
   5.00 8911 1780
   5.57 2701 484
   6.18 6533 1056
   6.20 5461 880
   6.41 1541 240
   8.22 8321 1012
   8.33 4033 484
   8.52 2047 240
   8.65 1387 160

Before n is strong pseduoprime, it must be a Fermat pseduoprime.
Let p1, ... pk be distinct prime factors of N:

Number of bases (mod N) = gcd(p1-1, N-1) ... gcd(pk-1, N-1)

So, picking top 2 numbers from the list:
1891 = 31 * 61, gcd(30, 1890) * gcd(60, 1890) = 30 * 30 = 900
8911 = 7*19*67, gcd(6, 8910) * gcd(18, 8910) * gcd(66, 8910) = 6 * 18 * 66 = 7128

No wonder so many strong pseudoprimes.

Gerald H's example, N = 803189 * 485909
gcd(N-1, 803188) * gcd(N-1, 485908) = 4*4 = 16
Removing 2 trivial bases of ±1, you are down to 14
Some Fermat pseudoprimes might not survive the strong test, so possibly few strong pseudoprimes.
Find all posts by this user
Quote this message in a reply
02-17-2019, 04:15 AM
Post: #40
RE: [HP35s] Program for prime number (brute force)
Meanwhile I extended the search a bit:
Code:
       r     n       k
       12.02 1563151 130048
       12.02 1844267 153456
       16.02 1911397 119284
       16.03 1600117 99844
       30.05 1968557 65520
       32.05 1663213 51892
       32.06 1389581 43348
       35.77 1922801 53748
       40.06 2078737 51892
       40.07 1333603 33280
       43.71 1777793 40676
       53.41 2102437 39364
       60.09 1999759 33280
       60.10 1589531 26448
       60.11 1275347 21216
       70.10 2121013 30256
       70.12 1546021 22048
       84.15 1455451 17296
       93.50 1459009 15604
       …

For \(1563151=1021×1531\) or \(1844267=1109×1663\) it's realistic to hit a strong liar after a few attempts.

Cheers
Thomas
Find all posts by this user
Quote this message in a reply
Post Reply 




User(s) browsing this thread: