\ FILE : pipi.frt
\ LANGUAGE : ANSI Forth
\ COPYRIGHT : Albert van der Horst FIG Chapter Holland
\ DATE : 1996-03-14
\ LAST CHANGE : 1996-03-14
\ This program may be copied and distributed freely as is.
\ It may be modified only for
\ - algorithmic improvements
\ - porting to other Forth systems.
\ Expressly forbidden is
\ 1. gratitious de-ansification, like use of non-standard comment symbols
\ 2. lower casing
\ 3. comment stripping
\ 4. addition of system-specific words like source control tools such as
\ MARKER look-alikes
\ ?Prime tests whether the single precision number P is prime
\ Cases 0 1 are in fact illegal but return TRUE
: ?Prime ( n P -- flag )
LOCAL P
P 4 U< IF TRUE EXIT THEN \ Prevent silly infinite loop
P 1 AND 0= IF FALSE EXIT THEN \ Handle even numbers other than 2
P 3 DO
P I /MOD
I < IF DROP TRUE LEAVE THEN
0= IF FALSE LEAVE THEN
2 +LOOP
;
\ N3 are the amount of numbers <= N1 that are dismissed by the prime N2,
\ i.e. it is divisible by N2 but by no smaller prime.
\ Requires N2<=N1
: DISMISS ( d N1, n N2 -- d N3 )
\ 2DUP CR ." Dismissing " . .
LOCAL P DLOCAL N
N 1 P M*/ DLOCAL N'
N' P S>D D< IF 1. EXIT THEN \ Only P itself
N'
P 2 ?DO
I ?Prime IF
N' I RECURSE D- \ Dismissed by a smaller prime
THEN
LOOP
;
\ Return PI(N1) i.e. the number of primes <= N1
: PI ( d N1 -- d N2 )
DLOCAL N
N \ All numbers are prime candidates
1. D- \ Exclude 1
N D2/ D- \ Exclude multiples of 2
1. D+ \ Except 2 itself
0 3 DO \ Upper limit never reached
I I UM* N D> IF LEAVE THEN
I ?Prime IF
N I DISMISS D- \ Exclude multiples of I
1. D+ \ Except I itself
THEN
2 +LOOP
;
: findproblem
1. DLOCAL OLD
0 3 DO
I S>D PI 2DUP OLD D= 0= \ Is I prime according to PI?
I ?Prime 0= 0= \ Is I really a prime?
<> IF DROP ." Wrong: " I . LEAVE THEN
TO OLD
KEY? IF ." tested till " I . LEAVE THEN
LOOP
;