File Coverage

blib/lib/Math/NumSeq/Powerful.pm

Criterion	Covered	Total	%
statement	49	61	80.3
branch	12	24	50.0
condition	3	9	33.3
subroutine	10	12	83.3
pod	4	4	100.0
total	78	110	70.9

# Copyright 2011, 2012, 2013, 2014, 2015 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 .

package Math::NumSeq::Powerful;
use 5.004;
use strict;
use List::Util 'min';

use vars '$VERSION','@ISA';
$VERSION = 72;
use Math::NumSeq 7; # v.7 for _is_infinite()
use Math::NumSeq::Base::IteratePred;
@ISA = ('Math::NumSeq::Base::IteratePred',
        'Math::NumSeq');
*_is_infinite = \&Math::NumSeq::_is_infinite;

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


# use constant name => Math::NumSeq::__('Powerful Numbers');
use constant i_start => 1;

sub description {
  my ($self) = @_;
  if (ref $self) {
    return ("Integers with $self->{'powerful_type'} prime factors "
            . ($self->{'power'} == 2 ? "squared"
               : $self->{'power'} == 3 ? "cubed"
               : "$self->{'power'}th power")
            . " or higher.");
  } else {
    return Math::NumSeq::__('Integers which some prime factors squared or higher.');
  }
}

use constant parameter_info_array =>
  [
   { name    => 'powerful_type',
     type    => 'enum',
     default => 'some',
     choices => ['some','all'],
     choices_display => [Math::NumSeq::__('Some'),
                         Math::NumSeq::__('All'),
                        ],
     # description => Math::NumSeq::__(''),
   },
   { name    => 'power',
     type    => 'integer',
     default => '2',
     minimum => 2,
     width   => 2,
     # description => Math::NumSeq::__(''),
   },
  ];

sub values_min {
  my ($self) = @_;
  return ($self->{'powerful_type'} eq 'some'
          ? 2 ** $self->{'power'}
          : 1);
}

# cf A168363 squares and cubes of primes, the primitives of "all" power=2
#    A112526 0/1 charactistic of A001694 all k>=2
#    A005934 highly powerful - new highest product of exponents
#    A005117 the squarefrees, all k<2, complement of some k>=2
#    A001597 perfect powers m^k for k>=2
#    A052485 weak, non-powerful, some k<2
#
my %oeis_anum = (some => [undef,
                          undef,
                          'A013929', # 2 non squarefree, divisible by some p^2
                          'A046099', # 3 non cubefree, divisible by some p^3
                          'A046101', # 4 non 4th-free
                          # OEIS-Catalogue: A013929
                          # OEIS-Catalogue: A046099 power=3
                          # OEIS-Catalogue: A046101 power=4
                         ],
                 all  => [undef,
                          undef,
                          'A001694', # 2 all p^k has k >= 2
                          'A036966', # 3 all p^k has k >= 3
                          'A036967', # 4 all p^k has k >= 4
                          'A069492', # 5 all p^k has k >= 5
                          'A069493', # 6
                          # OEIS-Catalogue: A001694 powerful_type=all
                          # OEIS-Catalogue: A036966 powerful_type=all power=3
                          # OEIS-Catalogue: A036967 powerful_type=all power=4
                          # OEIS-Catalogue: A069492 powerful_type=all power=5
                          # OEIS-Catalogue: A069493 powerful_type=all power=6
                         ],
                );
sub oeis_anum {
  my ($self) = @_;
  return $oeis_anum{$self->{'powerful_type'}}->[$self->{'power'}];
}

sub pred {
  my ($self, $value) = @_;
  ### SquareFree pred(): $value

  $value = abs($value);
  unless ($value >= 0) {
    return 0;
  }
  if ($value <= 0xFFFF_FFFF) {
    $value = "$value"; # numize Math::BigInt for speed
  }

  if ($value < 1 || $value != int($value)) {
    return 0;
  }

  my $power = $self->{'power'};
  my $limit = "$value" ** (1/$power) + 3;
  my $reduced_limit = min($limit,65535);

  my $p = 2;
  for ( ; $p <= $reduced_limit; $p += 2-($p==2)) {
    next if ($value % $p);
    ### prime factor: $p

    $value /= $p;
    my $count = 1;
    while (($value % $p) == 0) {
      ++$count;
      if ($count >= $power && $self->{'powerful_type'} eq 'some') {
        ### found some factor of desired power ...
        return 1;
      }
      $value /= $p;
    }
    if ($count < $power && $self->{'powerful_type'} ne 'some') {
      ### all/all_prim prime without desired power ...
      return 0;
    }

    $limit = "$value" ** (1/$power) + 3;
    $reduced_limit = min($limit,65535);
    ### divided out: "$p, now value=$value, new limit $limit reduced $reduced_limit"
  }
  ### final value: $value

  if ($p < $limit) {
    ### value too big to check ...
    return undef;
  }
  if ($self->{'powerful_type'} eq 'some') {
    ### some, no suitable power found ...
    return 0;
  } else {
    ### all, ok if reduced to value==1 ...
    return ($value == 1);
  }
}

1;
__END__