# Copyright 2012, 2013, 2014 Kevin Ryde
# This file is part of Math-NumSeq.
#
# Math-NumSeq is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by the
# Free Software Foundation; either version 3, or (at your option) any later
# version.
#
# Math-NumSeq is distributed in the hope that it will be useful, but
# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
# or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
# for more details.
#
# You should have received a copy of the GNU General Public License along
# with Math-NumSeq. If not, see <http://www.gnu.org/licenses/>.
package Math::NumSeq::EratosthenesStage;
use 5.004;
use strict;
use vars '$VERSION', '@ISA';
$VERSION = 72;
use Math::NumSeq;
use List::Util 'max';
@ISA = ('Math::NumSeq');
*_is_infinite = \&Math::NumSeq::_is_infinite;
use Math::NumSeq::Primes;
# uncomment this to run the ### lines
#use Smart::Comments;
# use constant name => Math::NumSeq::__('...');
use constant description => Math::NumSeq::__('Stages of the sieve of Eratosthenes.');
use constant i_start => 1;
use constant characteristic_smaller => 0;
use constant characteristic_integer => 1;
use constant parameter_info_array =>
[
{
name => 'stage',
type => 'integer',
default => '3',
width => 4,
# description => Math::NumSeq::__('...'),
},
];
use constant values_min => 2;
#------------------------------------------------------------------------------
# cf A179546 first a(n) numbers killed
# A179545 sum jump
# A066680 sieve up to p^2 only
# A083140 diagonals of stages
my @oeis_anum = (undef, # 0, all integers 2 upwards
'A004280', # 1
'A038179', # 2
);
sub oeis_anum {
my ($self) = @_;
return $oeis_anum[$self->{'stage'}];
}
#------------------------------------------------------------------------------
sub new {
my $self = shift->SUPER::new(@_);
my $stage = $self->{'stage'};
my @primes;
my $seq = Math::NumSeq::Primes->new;
$seq->next; # not 2
foreach (2 .. $stage) {
my ($i, $value) = $seq->next;
push @primes, $value;
}
$self->{'primes'} = \@primes;
### @primes
return $self;
}
sub rewind {
my ($self) = @_;
$self->{'i'} = $self->i_start;
$self->{'a'} = 2;
$self->{'size'} = 5;
$self->{'lo'} = 0;
$self->{'hi'} = 0;
}
sub next {
my ($self) = @_;
my $flags = $self->{'flags'};
my $lo = $self->{'lo'};
my $hi = $self->{'hi'};
for (;;) {
my $a = $self->{'a'}++;
### $a
if ($a > $hi) {
my $size = $self->{'size'} = int (1.08 * $self->{'size'});
# lo to hi inclusive,
# flags 0 to size=hi-lo inclusive
$lo = $self->{'lo'} = $a;
$hi = $self->{'hi'} = $lo + $size;
my $flags_str = '';
vec($flags_str,$size,1) = 0; # pre-extend
$flags = $self->{'flags'} = \$flags_str;
if ($self->{'stage'} > 0) {
for (my $n = max (4-$lo, -$lo % 2); $n <= $size; $n += 2) {
### flag even: "n=$n is ".($n+$lo)
vec($flags_str,$n,1) = 1;
}
foreach my $prime (@{$self->{'primes'}}) {
my $p2 = 2*$prime;
for (my $n = max (3*$prime - $lo,
-$lo % $prime);
$n <= $size;
$n += $p2) {
vec($flags_str,$n,1) = 1;
}
}
}
}
if (vec($$flags,$a-$lo,1)) {
### skip: $a
next;
}
### return: "i=$self->{'i'} a=$a"
return ($self->{'i'}++, $a);
}
}
1;
__END__
stage 0: 2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17
stage 1: 2,3,5,7,9,11,13,15,17
stage 2: 2,3,5,7,11,13,17