! Two routines related to finding prime numbers.
! Keywords: Prime number;
PROCEDURE PRIMETRY:
! Test whether PRMCAND8 is a prime number, i.e.,
not divisible by any other integer except 1 and PRMCAND8;
! Inputs:
PRMCAND8 = an integer to be tested if it is a prime number.
Outputs:
PRMFACTOR = largest factor of PRMCAND8
= 1 if PRMCAND8 is a prime,
= 0 if cannot tell, e.g., PRMCAND8 is too large,
= -1 if not an integer;
PRMABS = @ABS( PRMCAND8);
@IFC( PRMABS #GT# 2^46:
PRMFACTOR = 0; ! Too big to test;
@ELSE
@IFC( @FLOOR( PRMABS) #NE# PRMABS:
PRMFACTOR = -1; ! Not an integer?;
@ELSE
PRMFACTOR = 1; ! Assume prime until proven otherwise;
@IFC( PRMABS #NE# 2 #AND# PRMABS #NE# 3:
! Test if 2 or 3 is a factor of PRMCAND8;
@IFC( PRMABS #EQ# 2*@FLOOR(PRMABS/2):
PRMFACTOR = 2;
@ELSE
@IFC( PRMABS #EQ# 3*@FLOOR(PRMABS/3):
PRMFACTOR = 3;
);
);
! Test if 5, 7, 11, 13, 17, 19, ... are factors;
PRMDIVISOR = 5;
PRMSKIP = 2;
@WHILE( PRMDIVISOR*PRMDIVISOR #LE# PRMABS #AND# PRMFACTOR #EQ# 1:
@IFC( PRMDIVISOR*@FLOOR(PRMABS/PRMDIVISOR) #EQ# PRMABS:
PRMFACTOR = PRMDIVISOR;
);
PRMDIVISOR = PRMDIVISOR + PRMSKIP;
PRMSKIP = 6 - PRMSKIP;
);
);
);
);
ENDPROCEDURE
PROCEDURE PRIMELARGEST:
! Find largest prime number <= PRMUB;
! Inputs:
PRMUB = upper bound on prime we are searching for,
Outputs:
PRMLARGE = largest prime <= PRMUB
;
PRMFACTOR = 2; ! Not prime/has a factor until we find a prime;
PRMCAND8 = @ABS(PRMUB); ! Assume meant positive integers;
! Cannot achieve perfect precision for large numbers, so restrict size;
PRMCAND8 = @SMIN( 2^46, PRMCAND8);
PRMCAND8 = @FLOOR( PRMCAND8); ! Consider only integers;
@IFC( PRMCAND8 #GT# 2: ! 2 is special;
! Consider only odd integers;
@IFC( 2*@FLOOR( PRMCAND8/2) #EQ# PRMCAND8: PRMCAND8 = PRMCAND8 - 1);
PRMCAND8 = PRMCAND8 + 2;
@WHILE( PRMFACTOR #GT# 1:
PRMCAND8 = PRMCAND8 - 2; ! Next candidate to test for primality;
! @write( ' Test primality of ', PRMCAND8, @newline(1));
PRIMETRY; ! Test PRMCAND8 for primality;
);
);
PRMLARGE = PRMCAND8;
ENDPROCEDURE
CALC:
PRMUB = 84893256347; ! A number to test for primality;
PRIMELARGEST; ! Test for primality;
@WRITE( PRMLARGE, ' = largest prime <= ', PRMUB, @NEWLINE(1));
ENDCALC
|