src/wheel-gen.pl

   1 #!/usr/bin/perl -w
   2
   3 eval 'exec /usr/bin/perl -S $0 ${1+"$@"}'
   4   if 0;
   5
   6 use strict;
   7 (my $ME = $0) =~ s|.*/||;
   8
   9 # A global destructor to close standard output with error checking.
  10 sub END
  11 {
  12   defined fileno STDOUT
  13     or return;
  14   close STDOUT
  15     and return;
  16   warn "$ME: closing standard output: $!\n";
  17   $? ||= 1;
  18 }
  19
  20 sub is_prime ($)
  21 {
  22   my ($n) = @_;
  23   use integer;
  24
  25   $n == 2
  26     and return 1;
  27
  28   my $d = 2;
  29   my $w = 1;
  30   while (1)
  31     {
  32       my $q = $n / $d;
  33       $n == $q * $d
  34         and return 0;
  35       $d += $w;
  36       $q < $d
  37         and last;
  38       $w = 2;
  39     }
  40   return 1;
  41 }
  42
  43 {
  44   @ARGV == 1
  45     or die "$ME: missing argument\n";
  46
  47   my $wheel_size = $ARGV[0];
  48
  49   my @primes = (2);
  50   my $product = $primes[0];
  51   my $n_primes = 1;
  52   for (my $i = 3; ; $i += 2)
  53     {
  54       if (is_prime $i)
  55         {
  56           push @primes, $i;
  57           $product *= $i;
  58           ++$n_primes == $wheel_size
  59             and last;
  60         }
  61     }
  62
  63   my $ws_m1 = $wheel_size - 1;
  64   print <<EOF;
  65 /* The first $ws_m1 elements correspond to the incremental offsets of the
  66    first $wheel_size primes (@primes).  The $wheel_size(th) element is the
  67    difference between that last prime and the next largest integer
  68    that is not a multiple of those primes.  The remaining numbers
  69    define the wheel.  For more information, see
  70    http://www.utm.edu/research/primes/glossary/WheelFactorization.html.  */
  71 EOF
  72
  73   my @increments;
  74   my $prev = 2;
  75   for (my $i = 3; ; $i += 2)
  76     {
  77       my $rel_prime = 1;
  78       foreach my $divisor (@primes)
  79         {
  80           $i != $divisor && $i % $divisor == 0
  81             and $rel_prime = 0;
  82         }
  83
  84       if ($rel_prime)
  85         {
  86           #warn $i, ' ', $i - $prev, "\n";
  87           push @increments, $i - $prev;
  88           $prev = $i;
  89
  90           $product + 1 < $i
  91             and last;
  92         }
  93     }
  94
  95   print join (",\n", @increments), "\n";
  96
  97   exit 0;
  98 }