File Coverage

blib/lib/Math/NumSeq/Powerful.pm
Criterion Covered Total %
statement 50 62 80.6
branch 12 24 50.0
condition 3 9 33.3
subroutine 10 12 83.3
pod 4 4 100.0
total 79 111 71.1


line stmt bran cond sub pod time code
1             # Copyright 2011, 2012, 2013, 2014 Kevin Ryde
2              
3             # This file is part of Math-NumSeq.
4             #
5             # Math-NumSeq is free software; you can redistribute it and/or modify
6             # it under the terms of the GNU General Public License as published by the
7             # Free Software Foundation; either version 3, or (at your option) any later
8             # version.
9             #
10             # Math-NumSeq is distributed in the hope that it will be useful, but
11             # WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
12             # or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
13             # for more details.
14             #
15             # You should have received a copy of the GNU General Public License along
16             # with Math-NumSeq. If not, see .
17              
18             package Math::NumSeq::Powerful;
19 2     2   34980 use 5.004;
  2         9  
  2         80  
20 2     2   11 use strict;
  2         3  
  2         139  
21 2     2   12 use List::Util 'min';
  2         4  
  2         201  
22              
23 2     2   10 use vars '$VERSION','@ISA';
  2         4  
  2         118  
24             $VERSION = 71;
25 2     2   523 use Math::NumSeq 7; # v.7 for _is_infinite()
  2         38  
  2         44  
26 2     2   1185 use Math::NumSeq::Base::IteratePred;
  2         4  
  2         160  
27             @ISA = ('Math::NumSeq::Base::IteratePred',
28             'Math::NumSeq');
29             *_is_infinite = \&Math::NumSeq::_is_infinite;
30              
31             # uncomment this to run the ### lines
32             # use Smart::Comments;
33              
34              
35             # use constant name => Math::NumSeq::__('Powerful Numbers');
36 2     2   11 use constant i_start => 1;
  2         4  
  2         380  
37              
38             sub description {
39 0     0 1 0 my ($self) = @_;
40 0 0       0 if (ref $self) {
41 0 0       0 return ("Integers with $self->{'powerful_type'} prime factors "
    0          
42             . ($self->{'power'} == 2 ? "squared"
43             : $self->{'power'} == 3 ? "cubed"
44             : "$self->{'power'}th power")
45             . " or higher.");
46             } else {
47 0         0 return Math::NumSeq::__('Integers which some prime factors squared or higher.');
48             }
49             }
50              
51 2         13 use constant parameter_info_array =>
52             [
53             { name => 'powerful_type',
54             type => 'enum',
55             default => 'some',
56             choices => ['some','all'],
57             choices_display => [Math::NumSeq::__('Some'),
58             Math::NumSeq::__('All'),
59             ],
60             # description => Math::NumSeq::__(''),
61             },
62             { name => 'power',
63             type => 'integer',
64             default => '2',
65             minimum => 2,
66             width => 2,
67             # description => Math::NumSeq::__(''),
68             },
69 2     2   10 ];
  2         4  
70              
71             sub values_min {
72 3     3 1 8 my ($self) = @_;
73 3 100       20 return ($self->{'powerful_type'} eq 'some'
74             ? 2 ** $self->{'power'}
75             : 1);
76             }
77              
78             # cf A168363 squares and cubes of primes, the primitives of "all" power=2
79             # A112526 0/1 charactistic of A001694 all k>=2
80             # A005934 highly powerful - new highest product of exponents
81             # A005117 the squarefrees, all k<2, complement of some k>=2
82             # A001597 perfect powers m^k for k>=2
83             # A052485 weak, non-powerful, some k<2
84             #
85             my %oeis_anum = (some => [undef,
86             undef,
87             'A013929', # 2 non squarefree, divisible by some p^2
88             'A046099', # 3 non cubefree, divisible by some p^3
89             'A046101', # 4 non 4th-free
90             # OEIS-Catalogue: A013929
91             # OEIS-Catalogue: A046099 power=3
92             # OEIS-Catalogue: A046101 power=4
93             ],
94             all => [undef,
95             undef,
96             'A001694', # 2 all p^k has k >= 2
97             'A036966', # 3 all p^k has k >= 3
98             'A036967', # 4 all p^k has k >= 4
99             'A069492', # 5 all p^k has k >= 5
100             'A069493', # 6
101             # OEIS-Catalogue: A001694 powerful_type=all
102             # OEIS-Catalogue: A036966 powerful_type=all power=3
103             # OEIS-Catalogue: A036967 powerful_type=all power=4
104             # OEIS-Catalogue: A069492 powerful_type=all power=5
105             # OEIS-Catalogue: A069493 powerful_type=all power=6
106             ],
107             );
108             sub oeis_anum {
109 0     0 1 0 my ($self) = @_;
110 0         0 return $oeis_anum{$self->{'powerful_type'}}->[$self->{'power'}];
111             }
112              
113             sub pred {
114 1002     1002 1 40743 my ($self, $value) = @_;
115             ### SquareFree pred(): $value
116              
117 1002         1028 $value = abs($value);
118 1002 50       1819 unless ($value >= 0) {
119 0         0 return 0;
120             }
121 1002 100       1974 if ($value <= 0xFFFF_FFFF) {
122 1000         1156 $value = "$value"; # numize Math::BigInt for speed
123             }
124              
125 1002 50 33     4484 if ($value < 1 || $value != int($value)) {
126 0         0 return 0;
127             }
128              
129 1002         1607 my $power = $self->{'power'};
130 1002         2751 my $limit = "$value" ** (1/$power) + 3;
131 1002         1874 my $reduced_limit = min($limit,65535);
132              
133 1002         1254 my $p = 2;
134 1002         1986 for ( ; $p <= $reduced_limit; $p += 2-($p==2)) {
135 2031 100       5638 next if ($value % $p);
136             ### prime factor: $p
137              
138 424         1036 $value /= $p;
139 424         1194 my $count = 1;
140 424         792 while (($value % $p) == 0) {
141 1306         50451 ++$count;
142 1306 50 33     4877 if ($count >= $power && $self->{'powerful_type'} eq 'some') {
143             ### found some factor of desired power ...
144 0         0 return 1;
145             }
146 1306         2578 $value /= $p;
147             }
148 424 50 33     2013 if ($count < $power && $self->{'powerful_type'} ne 'some') {
149             ### all/all_prim prime without desired power ...
150 0         0 return 0;
151             }
152              
153 424         992 $limit = "$value" ** (1/$power) + 3;
154 424         1551 $reduced_limit = min($limit,65535);
155             ### divided out: "$p, now value=$value, new limit $limit reduced $reduced_limit"
156             }
157             ### final value: $value
158              
159 1002 50       1870 if ($p < $limit) {
160             ### value too big to check ...
161 0         0 return undef;
162             }
163 1002 50       2113 if ($self->{'powerful_type'} eq 'some') {
164             ### some, no suitable power found ...
165 0         0 return 0;
166             } else {
167             ### all, ok if reduced to value==1 ...
168 1002         3765 return ($value == 1);
169             }
170             }
171              
172             1;
173             __END__