1#!/usr/bin/perl -w 2# Generate the spokes of a wheel, for wheel factorization. 3 4# Copyright (C) 2001, 2005, 2009-2010 Free Software Foundation, Inc. 5 6# This program is free software: you can redistribute it and/or modify 7# it under the terms of the GNU General Public License as published by 8# the Free Software Foundation, either version 3 of the License, or 9# (at your option) any later version. 10 11# This program is distributed in the hope that it will be useful, 12# but WITHOUT ANY WARRANTY; without even the implied warranty of 13# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 14# GNU General Public License for more details. 15 16# You should have received a copy of the GNU General Public License 17# along with this program. If not, see <http://www.gnu.org/licenses/>. 18 19eval 'exec /usr/bin/perl -S $0 ${1+"$@"}' 20 if 0; 21 22use strict; 23(my $ME = $0) =~ s|.*/||; 24 25# A global destructor to close standard output with error checking. 26sub END 27{ 28 defined fileno STDOUT 29 or return; 30 close STDOUT 31 and return; 32 warn "$ME: closing standard output: $!\n"; 33 $? ||= 1; 34} 35 36sub is_prime ($) 37{ 38 my ($n) = @_; 39 use integer; 40 41 $n == 2 42 and return 1; 43 44 my $d = 2; 45 my $w = 1; 46 while (1) 47 { 48 my $q = $n / $d; 49 $n == $q * $d 50 and return 0; 51 $d += $w; 52 $q < $d 53 and last; 54 $w = 2; 55 } 56 return 1; 57} 58 59{ 60 @ARGV == 1 61 or die "$ME: missing argument\n"; 62 63 my $wheel_size = $ARGV[0]; 64 65 my @primes = (2); 66 my $product = $primes[0]; 67 my $n_primes = 1; 68 for (my $i = 3; ; $i += 2) 69 { 70 if (is_prime $i) 71 { 72 push @primes, $i; 73 $product *= $i; 74 ++$n_primes == $wheel_size 75 and last; 76 } 77 } 78 79 my $ws_m1 = $wheel_size - 1; 80 print <<EOF; 81/* The first $ws_m1 elements correspond to the incremental offsets of the 82 first $wheel_size primes (@primes). The $wheel_size(th) element is the 83 difference between that last prime and the next largest integer 84 that is not a multiple of those primes. The remaining numbers 85 define the wheel. For more information, see 86 http://www.utm.edu/research/primes/glossary/WheelFactorization.html. */ 87EOF 88 89 my @increments; 90 my $prev = 2; 91 for (my $i = 3; ; $i += 2) 92 { 93 my $rel_prime = 1; 94 foreach my $divisor (@primes) 95 { 96 $i != $divisor && $i % $divisor == 0 97 and $rel_prime = 0; 98 } 99 100 if ($rel_prime) 101 { 102 #warn $i, ' ', $i - $prev, "\n"; 103 push @increments, $i - $prev; 104 $prev = $i; 105 106 $product + 1 < $i 107 and last; 108 } 109 } 110 111 print join (",\n", @increments), "\n"; 112 113 exit 0; 114} 115