File Coverage

blib/lib/Math/NumSeq/WoodallNumbers.pm
Criterion Covered Total %
statement 89 96 92.7
branch 21 26 80.7
condition 7 9 77.7
subroutine 18 19 94.7
pod 7 7 100.0
total 142 157 90.4


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::WoodallNumbers;
19 2     2   6219 use 5.004;
  2         5  
20 2     2   8 use strict;
  2         6  
  2         35  
21              
22 2     2   7 use vars '$VERSION', '@ISA';
  2         3  
  2         88  
23             $VERSION = 72;
24              
25 2     2   352 use Math::NumSeq;
  2         3  
  2         67  
26             @ISA = ('Math::NumSeq');
27             *_is_infinite = \&Math::NumSeq::_is_infinite;
28              
29 2     2   340 use Math::NumSeq::Fibonacci;
  2         2  
  2         70  
30             *_blog2_estimate = \&Math::NumSeq::Fibonacci::_blog2_estimate;
31              
32             # uncomment this to run the ### lines
33             #use Smart::Comments;
34              
35              
36             # use constant name => Math::NumSeq::__('Woodall Numbers');
37 2     2   7 use constant description => Math::NumSeq::__('Woodall numbers i*2^i-1.');
  2         3  
  2         4  
38 2     2   6 use constant characteristic_increasing => 1;
  2         3  
  2         76  
39 2     2   6 use constant characteristic_integer => 1;
  2         2  
  2         66  
40 2     2   6 use constant values_min => 1;
  2         2  
  2         64  
41 2     2   7 use constant i_start => 1; # from 1*2^1-1==1
  2         1  
  2         70  
42              
43             # cf A002234 - n of Woodall primes
44             # A050918 - Woodall primes values
45             # A056821 - totient(woodall)
46             #
47 2     2   7 use constant oeis_anum => 'A003261';
  2         1  
  2         1087  
48              
49             # pow*i+1
50             my $uv_i_limit = do {
51             my $max = ~0;
52             my $limit = 1;
53              
54             for (my $i = 1; $i < 1000; $i++) {
55             if ($i <= (($max>>1)+1) >> ($i-1)) {
56             $limit = $i;
57             } else {
58             last;
59             }
60             }
61              
62             ### max : $max
63             ### woodall: (1<<$limit)*$limit+1
64             ### assert: $limit*(1<<$limit)+1 <= $max
65              
66             $limit
67             };
68             ### $uv_i_limit
69              
70             sub rewind {
71 4     4 1 3 my ($self) = @_;
72 4         10 my $i = $self->{'i'} = $self->i_start;
73 4         7 $self->{'power'} = 2 ** $i;
74             }
75             sub seek_to_i {
76 14     14 1 570 my ($self, $i) = @_;
77             ### seek_to_i(): $i
78 14         12 $self->{'i'} = $i;
79 14 100       22 if ($i >= $uv_i_limit) {
80 2         6 $i = Math::NumSeq::_to_bigint($i);
81             }
82 14         192 $self->{'power'} = 2 ** $i;
83             }
84             sub _UNTESTED__seek_to_value {
85 0     0   0 my ($self, $value) = @_;
86 0         0 $self->seek_to_i($self->value_to_i_ceil($value));
87             }
88              
89             # diff = (i+1)*2^(i+1)-1 - (i*2^i-1)
90             # = i*2^(i+1) + 2^(i+1) - i*2^i
91             # = i*(2^(i+1) - 2^i) + 2^(i+1)
92             # = i*2^i + 2^(i+1)
93             # 2*(i*2^i-1) + 2*2^i + 1
94             # = 2*i*2^i-2 + 2*2^i + 1
95             # = (2*i+2)*2^i - 1
96             # = (i+1)*2^(i+1) - 2
97             #
98             sub next {
99 25     25 1 376 my ($self) = @_;
100 25         30 my $i = $self->{'i'}++;
101 25 50       35 if ($i == $uv_i_limit) {
102 0         0 $self->{'power'} = Math::NumSeq::_to_bigint($self->{'power'});
103             }
104 25         31 my $value = $self->{'power'}*$i - 1;
105 25         415 $self->{'power'} *= 2;
106 25         172 return ($i, $value);
107             }
108              
109             sub ith {
110 14     14 1 390 my ($self, $i) = @_;
111 14         14 my $power;
112 14 100 66     42 if (! ref $i && $i >= $uv_i_limit) {
113 2         7 $power = Math::NumSeq::_to_bigint(1) << $i;
114             } else {
115 12         12 $power = 2**$i;
116             }
117 14         425 return $i * $power - 1;
118             }
119              
120             sub pred {
121 9     9 1 12 my ($self, $value) = @_;
122             ### WoodallNumbers pred(): $value
123              
124             {
125 9         5 my $int = int($value);
  9         9  
126 9 50       14 if ($value != $int) { return 0; }
  0         0  
127 9         8 $value = $int;
128             }
129 9 100 100     32 unless ($value >= 1
130             && ($value % 2) == 1) {
131 5         12 return 0;
132             }
133              
134 4         2 $value += 1; # now seeking $value == $exp * 2**$exp
135 4 50       10 if (_is_infinite($value)) {
136 0         0 return 0;
137             }
138              
139 4         5 my $exp = 0;
140 4         3 for (;;) {
141 12 100 66     27 if ($value <= $exp || $value % 2) {
142 4         11 return ($value == $exp);
143             }
144 8         7 $value >>= 1;
145 8         5 $exp++;
146             }
147             }
148              
149             # ENHANCE-ME: round_down_pow(value,2) then exp-=log2(exp) is maybe only
150             # 0,+1,+2 away from being correct
151             #
152             sub value_to_i_floor {
153 22526     22526 1 51656 my ($self, $value) = @_;
154             ### WoodallNumbers value_to_i_floor(): $value
155              
156 22526         13541 $value = int($value) + 1; # now seeking $value == $exp * 2**$exp
157 22526 100       23561 if ($value < 8) {
158 6         7 return 1;
159             }
160 22520 50       24806 if (_is_infinite($value)) {
161 0         0 return $value;
162             }
163              
164 22520         15745 my $i = 1;
165 22520         14145 my $power = ($value*0) + 2; # 2^1=2, inherit bignum 2
166              
167 22520         12346 for (;;) {
168 229351         128846 my $try_value = $i*$power;
169              
170             ### $i
171             ### try_value: "$try_value"
172              
173 229351 100       220563 if ($try_value >= $value) {
174 22520         26901 return $i - 1 + ($try_value == $value);
175             }
176              
177 206831 50       201912 if ($i == $uv_i_limit) {
178 0         0 $power = Math::NumSeq::_to_bigint($power);
179             }
180 206831         114597 $power *= 2;
181 206831         112885 $i++;
182             }
183             }
184              
185             # shared by Math::NumSeq::CullenNumbers
186             sub value_to_i_estimate {
187 19     19 1 248 my ($self, $value) = @_;
188              
189 19 100       22 if ($value < 1) {
190 8         9 return 0;
191             }
192              
193 11         78 my $i = _blog2_estimate($value);
194 11 100       235 if (! defined $i) {
195 10         10 $i = log($value) * (1/log(2));
196             }
197 11         13 $i -= log($i) * (1/log(2));
198 11         12 return int($i);
199             }
200              
201             1;
202             __END__