1.DATA STEP
data tmp(drop=i); n=2; output; do n=3 to 10000 by 2; do i=2 to n-1; if mod(n,i)=0 and i^=n-1 then leave; if i=n-1 then output; end; end; run; proc print;run;
2.macro from SAS_L
option mprint; %macro prime(n); data prime; do i=1 to &n; prime=1; do j=2 to ceil(i/2); if mod(i,j)=0 then prime=0; end; if prime=1 then output; end; run; %mend; %prime(10000);
%LET TOTAL = 2000; * Limits the number of prime numbers generated ; %LET DIM = 1000; * The size of the sieve arrays used ; %LET TIME = 30; * Time limit in seconds ; DATA _null_; ARRAY P{&DIM} _TEMPORARY_; * Prime numbers; ARRAY M{&DIM} _TEMPORARY_; * Multiples of prime numbers; TIMEOUT = DATETIME() + &TIME; * Time limit; FILE print NOTITLES; SQUARE = 4; DO X = 2 TO 10000; * Is X prime? ; IF DATETIME() >= TIMEOUT THEN STOP; * Time limit reached ; IF X = SQUARE THEN DO; * Extend sieve; IMAX + 1; IF IMAX >= &DIM THEN STOP; * Sieve size limit reached. ; SQUARE = M{IMAX + 1}; CONTINUE; END; * Find least prime factor (LPF). ; LPF = 0; DO I = 1 TO IMAX UNTIL (LPF); DO WHILE (M{I} < X); * Update sieve with new multiple. ; M{I} + P{I}; END; IF M{I} = X THEN LPF = P{I}; END; IF LPF THEN CONTINUE; * Composite number found. ; PUT X @; * Write prime number in output. ; N + 1; IF N >= &TOTAL THEN STOP; * Output maximum reached. ; ELSE IF N <= &DIM THEN DO; * Add prime number to sieve. ; P{N} = X; M{N} = X*X; END; END; STOP; RUN;
4.Rick Wicklin 's algorithm
/** The Sieve of Eratosthenes **/ proc iml; max = 10000; list = 2:max; /** find prime numbers in [2, max] **/ primes = j(1, ncol(list), .); /** allocate space to store primes **/ numPrimes = 0; p = list[1]; /** 2 is the first prime **/ do while (p##2 <= max); /** search through sqrt(max) **/ idx = loc( mod(list, p)=0 );/** find all numbers divisible by p **/ list = remove(list, idx); /** remove them. They are not prime. **/ numPrimes = numPrimes + 1; primes[numPrimes] = p; /** include p in the list of primes **/ p = list[1]; /** next prime number to sift **/ end; /** include remaining numbers greater than sqrt(max) **/ k = numPrimes; n = k + ncol(list); primes[k+1:n] = list; /** put them at the end of the list **/ primes = primes[ ,1:n]; /** get rid of excess storage space **/ print primes; quit;
5.discussions