File Coverage

blib/lib/Math/NumSeq/DeletablePrimes.pm

Criterion	Covered	Total	%
statement	70	70	100.0
branch	19	20	95.0
condition	6	6	100.0
subroutine	16	16	100.0
pod	4	4	100.0
total	115	116	99.1

# 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://mathworld.wolfram.com/DeletablePrime.html
#


package Math::NumSeq::DeletablePrimes;
use 5.004;
use strict;
use Math::Prime::XS 0.23 'is_prime'; # version 0.23 fix for 1928099

use vars '$VERSION', '@ISA';
$VERSION = 71;

use Math::NumSeq;
use Math::NumSeq::Primes;
@ISA = ('Math::NumSeq');
*_is_infinite = \&Math::NumSeq::_is_infinite;

# uncomment this to run the ### lines
#use Smart::Comments;


# use constant name => Math::NumSeq::__('Deletable Primes');
use constant description => Math::NumSeq::__('Deletable primes, being primes where deleting a digit gives another prime from which in turn a digit can be deleted, etc.');
use constant i_start => 1;
use constant characteristic_increasing => 1;
use constant values_min => 2;

use Math::NumSeq::Base::Digits
  'parameter_info_array';   # radix parameter


#------------------------------------------------------------------------------

# cf A080603 deleting a digit leaves a prime
#    A179336 primes containing at least one prime digit
#              a superset of the deletables
#
# http://mathworld.wolfram.com/TruncatablePrime.html
# A024770 left truncatable primes
# A024785 left truncatable primes, no zero digits
# A077390 left-and-right truncatable primes
# A137812 left or right truncatable primes
# 	 finite 149677 elements to 8939662423123592347173339993799
#
my @oeis_anum;

$oeis_anum[2] = 'A096246';
$oeis_anum[10] = 'A080608';
# OEIS-Catalogue: A096246 radix=2
# OEIS-Catalogue: A080608

sub oeis_anum {
  my ($self) = @_;
  return $oeis_anum[$self->{'radix'}];
}

#------------------------------------------------------------------------------


sub rewind {
  my ($self) = @_;
  $self->{'i'} = $self->i_start;
  my $primes_seq = $self->{'primes_seq'} = Math::NumSeq::Primes->new;
}

sub next {
  my ($self) = @_;
  my $primes_seq = $self->{'primes_seq'};
  for (;;) {
    (undef, my $prime) = $primes_seq->next
      or return;
    if (_prime_is_deletable($self,$prime)) {
      return ($self->{'i'}++, $prime);
    }
  }
}

sub pred {
  my ($self, $value) = @_;
  ### pred(): $value
  return ($self->Math::NumSeq::Primes::pred($value)
          && _prime_is_deletable($self,$value));
}

sub _prime_is_deletable {
  my ($self, @pending) = @_;

  my $radix = $self->{'radix'};
  my $target = ($radix == 2 ? 4 : $radix);

  while (@pending) {
    ### pending: join(', ',map{sprintf "%b", $_}@pending)

    my $value = pop @pending;
    next unless is_prime($value);
    if ($value < $target) {
      ### reached single digit, or 2,3 for binary ...
      return 1;
    }

    my $prev = -1;
    if ($radix == 10) {
      foreach my $i (($value =~ /^.(0+)/ ? length($1) : 0)
                     ..
                     length($value)-1) {
        my $digit = substr($value,$i,1);
        next if $digit eq $prev;
        $prev = $digit;
        push @pending, substr($value,0,$i) . substr($value,$i+1);
      }
    } else {

      for (my $pow = 1; $value > $pow; $pow *= $radix) {
        my $low = $value % $pow;
        my $high = int($value/$pow);

        my $digit = $high % $radix; # being deleted
        last if $high < $radix && $prev == 0; # don't delete high,0,...
        next if $digit eq $prev;
        $prev = $digit;


        push @pending, int($high/$radix) * $pow + $low;
      }
    }
  }
  return 0;
}

1;
__END__