Category: 
Fun Stuff 
Author/Contact Info 
thinker 
Description: 
A Sieve of Eratosthenes implemented with closures.
I was looking over an old Java project for Uni where I had implemented a Sieve as Objects connected by pipelines.
It occurred to me closures could do a similar job far more concisely. Here was my effort. I hope this is the correct section to post in. I found it fun, anyway. :)

#!/usr/bin/perl
# prints all prime numbers up to the value of the command line argumen
+t, or 500 if no args;
use strict;
use warnings;
use POSIX; # for ceil()
my @primes = ( 1, 2 ); # prime the list of primes :)
my $max = shift  500;
# initialise $max if no command line argument
my $enough = ceil ( $max ** 0.5 );
# we only need to filter until we have >= the sqrt
sub sieve {
# create a new closure to act as a part of the sieve
# the number that this closure will filter out multiples of
my $divisor = shift;
my $sieve_ref;
return sub {
my $number = shift;
# Does $number divide by divisor?
# If so, it's not prime, so bail out;
return unless $number % $divisor;
# if we have already created a next sieve,
# then let it deal with the number
# pass $number to it, and return;
if ( $sieve_ref ){
&$sieve_ref( $number );
return
};
#if we get this far, we have a new prime
# if we need to make a new sieve to filter,do so
$sieve_ref = sieve( $number ) if $number < $enough;
# push the prime onto the array
push @primes, $number;
};
};
# the main program
my $sieve = sieve( 2 ); # initialise
&$sieve( $_ ) for 3 .. $max;
print "The primes less than or equal to $max are @primes\n";

•Re: Sieve of Eratosthenes with closures
by merlyn (Sage) on Jul 20, 2003 at 23:07 UTC

Here's a slightly more direct implementation of that:
use strict;
$++;
my $generator = do {
my $num = 1;
sub { ++$num }
};
while (1) {
my $prime = $generator>();
print "$prime\n";
## don't let any multiple of these escape:
my $oldgenerator = $generator;
$generator = sub {
{
my $next_prime = $oldgenerator>();
redo if $next_prime % $prime == 0;
return $next_prime;
}
};
}
 Randal L. Schwartz, Perl hacker
Be sure to read my standard disclaimer if this is a reply.  [reply] [d/l] 
Re: Sieve of Eratosthenes with closures
by tilly (Archbishop) on Jul 21, 2003 at 00:14 UTC

Not all prime number sieves are the Sieve of Eratosthenes.
The key insight that makes the Sieve of Eratosthenes so nice is that for each prime you only have to do an operation on multiples of that prime. This causes the running time to be just slightly supralinear (IIRC it is O(n*log(log(n))) to get all the primes up to n) rather than being close to O(n**1.5) (yours should be O(n**1.5/log(n)).
Here is an actual Sieve of Eratosthenes implemented with closures. It isn't that fast, but it scales right. Memory usage should scale O(sqrt(n)).
#! /usr/bin/perl w
use strict;
sub build_sieve {
my $n = 0;
my @upcoming_factors = ();
my ($next_small_p, $sub_iter);
my $gives_next_square = 5;
return sub {
LOOP: {
if (not $n++) {
return 2; # Special case
}
if (not defined $upcoming_factors[0]) {
if ($n == $gives_next_square) {
if (not defined ($sub_iter)) {
# Be lazy to avoid an infinite loop...
$sub_iter = build_sieve();
$sub_iter>(); # Throw away 2
$next_small_p = $sub_iter>();
}
push @{$upcoming_factors[$next_small_p]}, $next_small_p;
$next_small_p = $sub_iter>();
my $next_p2 = $next_small_p * $next_small_p;
$gives_next_square = ($next_p2 + 1)/2;
shift @upcoming_factors;
redo LOOP;
}
else {
shift @upcoming_factors;
return 2*$n1;
}
}
else {
foreach my $i (@{$upcoming_factors[0]}) {
push @{$upcoming_factors[$i]}, $i;
}
shift @upcoming_factors;
redo LOOP;
}
}
}
}
# Produce as many primes as are asked for, or 100.
my $sieve = build_sieve();
print $sieve>(), "\n" for 1..(shift  100);
 [reply] [d/l] 

Wow, this looks like some complex stuff. How exactly would this be benneficial in a program? Or is it basically just concept and application?
 [reply] 

Closures allow you to organize code differently. The resulting style is as distinct from OO and procedural as they are from each other.
Right now Dominus is writing a book on the subject. You can find some information in this article on iterators, a number of my nodes talk about closures a bit, see Why I like functional programming for one of the better ones. Or you can search for more information.
Of course in this case the technique is just for fun. A somewhat shorter implementation of a sieve of Eratosthenes is:
sub sieve {
grep{@_[map$a*$_,$_..@_/($a=$_)]=0if$_[$_]>1}@_=0..pop
}
print join " ", sieve(100);
(Implementation due to tye at (tye)Re3: (Golf): Sieve of Eratosthenes.)  [reply] [d/l] 
Re: Sieve of Eratosthenes with closures
by djWestTexas (Initiate) on May 15, 2018 at 19:09 UTC

I know that this is a very old thread, but I just ran across it today while researching Sieve of Eratosthenes implementations. I am awed by the great programming, but much of it is beyond me. However, my main interest was to find a great SoE implementation for use and examples to others. I was really interested in outright performance. With this qualifier in mind, here are a few observations:
The fastest implementation was the c code by AnonymousMonk. However, it isn't perl and has memory limitations for finding larger primes. But it is blazing fast for small lists of primes  1st two million primes found in less than 1 second.
The fastest perl implementation was the code by tilly. It found all primes up to 2,000,003 in 1 second. It also scaled nearly linearly for up to 32 million primes (the last one I tested). Wow!
The original code by thinker was good. It found all primes up to 2,000,003 in 10 seconds. It didn't scale as well for larger values of primes.
The code by merlynn has some kind of problem. It took it 43 seconds to find all primes up to 2,000,003 and didn't scale well at all for larger values.
Hope my "years late" comment doesn't upset anyone's apple cart. But if anyone else, like me, runs across this thread, it may offer something useful.
If I need to find a bunch of primes quickly, I would use tilly's code!
djWestTexas  [reply] 

tilly's is a very nice closure implementation. It isn't the fastest pure Perl SoE though. For example, for the first 5761455 primes (primes up to 10^8) tilly's code takes 35.7s on my Macbook. The vector sieve on RosettaCode here takes 37.9s, and the string sieve here takes only 10.2 seconds, albeit with more memory and not an iterator. Using the purePerl segmented sieve in Math::Prime::Util::PP takes 7.1 seconds. All of these are terrifically slow compared to C.
If you're looking for even faster SoEs, I recommend Prime Sieve Benchmarks for a list of benchmarks. The C code from AnonymousMonk back in 2003 works fine, but is about 2x slower than the "Trivial 4line SoE" shown on that page, because it doesn't get the inner loop start correct. Of course they're all in C, though one is built into a Perl module. In general primesieve is the fastest and easiest solution for most projects.
To compare using XS modules, printing primes to 10^8. At this point a majority of our time is actually spent printing, which is why the specialized "print_primes()" routine stomps on everything else, since it has optimized C code (even C's printf is a major bottleneck once the sieving is fast enough).
 3.4s Math::Prime::XS
 2.8s Math::Prime::FastSieve
 2.0s Math::Prime::Util "say for @{primes(1e8)}"
 1.9s Math::Prime::Util "forprimes { say } 1e8"
 0.15s Math::Prime::Util "print_primes(1e8);"
Slightly better, here we're naively counting, except for the last one which uses a fast exact prime count method.
 1.7s Math::Prime::XS
 0.89s Math::Prime::FastSieve
 0.39s Math::Prime::Util "$n++ for @{primes(1e8)}"
 0.18s Math::Prime::Util "forprimes { $n++ } 1e8;"
 0.00026s Math::Prime::Util "say prime_count(1e8)"
 [reply] 

 [reply] 
Re: Sieve of Eratosthenes with closures
by hdb (Monsignor) on May 16, 2018 at 11:06 UTC

This is a great idea. I always wanted to implement a sieve that can just go on finding new primes, i.e. without specifying a maximum up to which the sieving will be conducted. This would require not only to store the primes up to a point but also kind of the largest multiple of each prime up to that point. Also, I think using the modulus operator % violates the sieving philosophy. I never was able to make it work as I got confused keeping the primes and the multiples.
Using closures is ideal for that. Every time one finds a new prime, one generates a new sifter that is responsible for that prime and to store the progress to date. Here is my formulation:
use strict;
use warnings;
sub sifter {
my $p = shift;
my $c = $p;
return sub { $c += $p while $c < $_[0]; return $c$_[0]; }
}
my @sieves;
my $n = 1;
loop: while( $n++ ){
$_>( $n ) or next loop for @sieves;
push @sieves, sifter( $n );
print "$n\n";
}
Clearly this is not a fast algorithm to find primes but I liked the idea...
 [reply] [d/l] [select] 
Re: Sieve of Eratosthenes with closures
by Anonymous Monk on Jul 21, 2003 at 12:51 UTC

Of course if you want fast, say 1 million primes in a second you would never use Perl. Perl is great for high level stuff. For raw algorithms.....
#include <stdio.h>
#include <malloc.h>
#include <time.h>
#define TEST(f,x) (*(f+(x)/16)&(1<<(((x)%16L)/2)))
#define SET(f,x) *(f+(x)/16)=1<<(((x)%16L)/2)
int
main(int argc, char *argv[])
{
unsigned char *ptr=NULL;
unsigned long max, mom, s=0, e=1;
register unsigned long test=1;
register unsigned long count=1;
time_t begin;
max = ( argc>1 ) ? atol(argv[1]) : 0xFFFFFFL;
ptr = (unsigned char *) calloc( ((max>>4)+1L), sizeof(char) );
if ( ! ptr ) {
printf( "Can't allocate enough memory!\n" );
system("pause");
return 1;
}
printf( "Searching prime numbers to : %ld\n", max );
begin = time (NULL);
while ( (test+=2) < max) if (!TEST(ptr, test) ) {
if ( ++count%0xFFL==0 ) {
printf ( "%ld prime number\x0d", count );
fflush(stdout);
}
for ( mom=3L*test; mom<max; mom+=test<<1 ) SET (ptr, mom);
}
printf( " %ld prime numbers found in %ld secs.\n\n", count, time(N
+ULL)begin);
printf( "Show prime numbers ( Enter end < start to exit )" );
while ( s<e ) {
printf ("\n\nStart of Area : ");
fflush (stdout); scanf ("%ld", &s);
printf ("End of Area : ");
fflush (stdout); scanf ("%ld", &e);
count = s2;
if ( s%2==0 ) count++;
while ((count+=2)<e) if (!TEST(ptr,count)) printf ("%ld\t", co
+unt);
}
free(ptr);
return 0;
}
/*
Searching prime numbers to : 16777215
1077871 prime numbers found in 1 secs.
Show prime numbers ( Enter end < start to exit )
Start of Area : 1
End of Area : 1000
1 3 5 7 11 13 17 19 23
+ 29
31 37 41 43 47 53 59 61 67
+ 71
73 79 83 89 97 101 103 107 109
+ 113
127 131 137 139 149 151 157 163 167
+ 173
179 181 191 193 197 199 211 223 227
+ 229
233 239 241 251 257 263 269 271 277
+ 281
283 293 307 311 313 317 331 337 347
+ 349
353 359 367 373 379 383 389 397 401
+ 409
419 421 431 433 439 443 449 457 461
+ 463
467 479 487 491 499 503 509 521 523
+ 541
547 557 563 569 571 577 587 593 599
+ 601
607 613 617 619 631 641 643 647 653
+ 659
661 673 677 683 691 701 709 719 727
+ 733
739 743 751 757 761 769 773 787 797
+ 809
811 821 823 827 829 839 853 857 859
+ 863
877 881 883 887 907 911 919 929 937
+ 941
947 953 967 971 977 983 991 997
Start of Area :
*/
 [reply] [d/l] 

