Talk:ALGOL 68-primes
Source code
Note the is probably prime PROC is based on code from The Miller–Rabin primality test task and the prelude.
<lang algol68># primes.incl.a68: prime related operators, procedure etc. #
# returns a sieve of primes up to n # OP PRIMESIEVE = ( INT n )[]BOOL: BEGIN [ 0 : n ]BOOL prime; prime[ 0 ] := prime[ 1 ] := FALSE; prime[ 2 ] := TRUE; FOR i FROM 3 BY 2 TO UPB prime DO prime[ i ] := TRUE OD; FOR i FROM 4 BY 2 TO UPB prime DO prime[ i ] := FALSE OD; FOR i FROM 3 BY 2 TO ENTIER sqrt( UPB prime ) DO IF prime[ i ] THEN FOR s FROM i * i BY i + i TO UPB prime DO prime[ s ] := FALSE OD FI OD; prime END; # PRIMESIEVE #
# modes and operators for extracting primes from a sieve #
INT first n primes extract op = 0; INT primes up to extract op = 1;
MODE PRIMEEXTRACT = STRUCT( INT extract op, INT n );
# returns a PRIMEEXTRACT with the first n primces extract op # OP EXTRACTFIRST = ( INT n )PRIMEEXTRACT: PRIMEEXTRACT( first n primes extract op, n );
# returns a PRIMEEXTRACT with primes up to extract op # OP EXTRACTPRIMESUPTO = ( INT n )PRIMEEXTRACT: PRIMEEXTRACT( primes up to extract op, n );
# returns an array of primes extracted from prime # # as specified by by the extract op and number of primes in ex # PRIO PRIMESFROMSIEVE = 9, FROMPRIMESIEVE = 9; OP PRIMESFROMSIEVE = ( PRIMEEXTRACT ex, []BOOL prime )[]INT: IF extract op OF ex = first n primes extract op THEN # extract the first n OF ex primes # # count the primes up n OF ex # INT p count := 0; FOR i TO UPB prime WHILE p count < n OF ex DO IF prime[ i ] THEN p count +:= 1 FI OD; # construct a list of the first n OF ex primes or p count, # # whichever is smaller # [ 1 : p count ]INT low prime; INT low pos := 0; FOR i WHILE low pos < UPB low prime DO IF prime[ i ] THEN low prime[ low pos +:= 1 ] := i FI OD; low prime ELSE # extract primes up to n OF ex primes # INT max prime = n OF ex; # count the primes up to max prime # INT p count := 0; FOR i TO max prime DO IF prime[ i ] THEN p count +:= 1 FI OD; # construct a list of the primes up to max prime # [ 1 : p count ]INT low prime; INT low pos := 0; FOR i WHILE low pos < UPB low prime DO IF prime[ i ] THEN low prime[ low pos +:= 1 ] := i FI OD; low prime FI # PRIMESFROMSIEVE # ; # alternative spelling of the operator name # OP FROMPRIMESIEVE = ( PRIMEEXTRACT ex, []BOOL prime )[]INT: ex PRIMESFROMSIEVE prime;
PROC is probably prime = ( LONG LONG INT n )BOOL: IF n = 2 THEN TRUE ELIF NOT ODD n THEN FALSE ELIF n < 3 THEN FALSE ELIF n < 10 000 THEN # smallish number = use trial division # BOOL is prime := TRUE; FOR i FROM 3 BY 2 TO ENTIER sqrt( SHORTEN SHORTEN n ) WHILE is prime := n MOD i /= 0 DO SKIP OD; is prime ELSE # miller rabin primality test # PROC pow mod = (LONG LONG INT b, in e, modulus)LONG LONG INT: BEGIN LONG LONG INT sq := b, e := in e; LONG LONG INT result := IF ODD e THEN b ELSE 1 FI; e OVERAB 2; WHILE e /= 0 DO sq := sq * sq MOD modulus; IF ODD e THEN result := result * sq MOD modulus FI ; e OVERAB 2 OD; result END # pow mod # ; INT k = 10; LONG LONG INT d := n - 1; INT s := 0; WHILE NOT ODD d DO d OVERAB 2; s +:= 1 OD; BOOL is prime := TRUE; TO k WHILE is prime DO LONG LONG INT a := 2 + ENTIER ( random * ( n - 3 ) ); LONG LONG INT x := pow mod( a, d, n ); IF x /= 1 THEN BOOL done := FALSE; TO s WHILE NOT done DO IF x = n - 1 THEN done := TRUE ELSE x := x * x MOD n FI OD; IF NOT done THEN IF x /= n-1 THEN is prime := FALSE FI FI FI OD; is prime FI # is probably prime # ;
- END primes.incl.a68 #</lang>