File Coverage

blib/lib/Math/Big/Factors.pm
Criterion Covered Total %
statement 97 100 97.0
branch 15 20 75.0
condition 5 10 50.0
subroutine 9 10 90.0
pod 2 2 100.0
total 128 142 90.1


line stmt bran cond sub pod time code
1             #############################################################################
2             # Math/Big/Factors.pm -- factor big numbers into prime factors
3              
4             package Math::Big::Factors;
5              
6             require 5.006002; # requires this Perl version or later
7              
8 2     2   58783 use strict;
  2         8  
  2         57  
9 2     2   10 use warnings;
  2         4  
  2         68  
10              
11 2     2   9 use Math::BigInt;
  2         3  
  2         12  
12 2     2   966 use Math::BigFloat;
  2         25127  
  2         12  
13 2     2   1101 use Math::Big;
  2         5  
  2         68  
14 2     2   10 use Exporter;
  2         4  
  2         1550  
15              
16             our $VERSION = '1.16';
17             our @ISA = qw( Exporter );
18             our @EXPORT_OK = qw( wheel factors_wheel
19             );
20              
21             sub wheel
22             {
23             # calculate a prime-wheel of order $o
24 12   50 12 1 1183 my $o = abs(shift || 1); # >= 1
25              
26             # some primitive wheels as shortcut:
27 12 100       31 return [ Math::BigInt->new(2), Math::BigInt->new(1) ] if $o == 1;
28              
29 9         37 my @primes = Math::Big::primes($o*5); # initial primes, get some more
30              
31 9         397 my $mul = Math::BigInt->new(1); my @wheel;
  9         383  
32 9         23 for (my $i = 0; $i < $o; $i++)
33             {
34             #print "$primes[$i]\n";
35 27         55 $mul *= $primes[$i]; push @wheel,$primes[$i];
  27         1264  
36             }
37             #print "Mul $mul\n";
38 9         14 my $last = $wheel[-1]; # get biggest initial
39             #print "last is $last\n";
40             # now sieve any number that is a multiply of one of the inital ones
41 9         33 @primes = (); # undef => leftover
42 9         15 foreach (@wheel)
43             {
44 27 100       1552 next if $_ == 2; # dont mark these, we skip 'em
45 18         1819 my $i = $_; my $add = $i;
  18         23  
46 18         38 while ($i < $mul)
47             {
48 462         30748 $primes[$i] = 1; $i += $add;
  462         8640  
49             }
50             }
51 9         711 push @wheel, Math::BigInt->new(1);
52 9         542 my $i = $last;
53 9         17 while ($i < $mul)
54             {
55 351 100       52258 push @wheel,$i if !defined $primes[$i]; $i += 2; # skip even ones
  351         6983  
56             }
57 9         1522 \@wheel;
58             }
59              
60             sub _transform_wheel
61             {
62             # from a given prime-wheel, calculate a increment table that can be used
63             # to step trough numbers
64             # input: ref to array with prime wheel
65             # output: ($restart,$ref_to_add_table);
66              
67 8     8   12 my (@wheel);
68 8         11 my $add = shift; shift @$add; # remove the first 2 from wheel
  8         12  
69              
70 8 100       27 if (@$add == 1) # order 1
71             {
72 2         4 my $two = Math::BigInt->new(2);
73             # (2,2) or (2,2,2,2,2,2) etc would do, too
74 2         87 @wheel = ($two->copy(),$two->copy(),$two->copy(),$two->copy());
75 2         113 return (1,\@wheel);
76             }
77             # from the list of divisors above create a add-table which we can take to
78             # increment from 3 onwards. The tabe consists of two parts, the second part
79             # will be repeatedly used
80 6         9 my $last = -1; my $mod = 2; my $i = 0;
  6         8  
  6         7  
81             # create the increment part for the initial primes (3,5, or 3,5,7 etc)
82 6         18 while ($add->[$i] != 1)
83             {
84 12         1187 $mod *= $add->[$i];
85 12 100       1142 push @wheel, $add->[$i] - $last if $last != -1; # skip the first
86             #print $wheel[-1],"\n" if $last != -1;
87 12         1181 $last = $add->[$i]; $i++;
  12         37  
88             }
89             #print "mod $mod\n";
90 6         585 my $border = $i-1; # account for ++
91 6         8 my $length = scalar @$add-$i;
92 6         10 my $ws = $border+$length; # remember this
93             #print "border: $border length $length $mod\n";
94              
95             # now we add two arrays in a row, both are equal except the first element
96             # which is in case A a step from the last inital prime to the second in list
97             # and in case B a step from '1' to the second in list
98              
99             #print "add[border+1]: ",$add->[$border+1]," add[border] $add->[$border]\n";
100 6         18 $wheel[$border] = $add->[$border+2]-$add->[$border];
101 6         585 $wheel[$border+$length] = $add->[$border+2]-1;
102             # and last add a wrap-around around $mod
103             #print "last: ",$add->[-1],"\n";
104 6         1029 $wheel[$border+$length-1] = 1+$mod-$add->[-1];
105 6         1380 $wheel[$border+$length*2-1] = $wheel[$border+$length-1];
106              
107 6         11 $i = $border + 1;
108             # now fill in the rest
109 6         17 while ($i < $length+$border-1)
110             {
111 104         201 $wheel[$i] = $add->[$i+2]-$add->[$i+1];
112 104         9389 $wheel[$i+$length] = $wheel[$i];
113 104         1092 $i++;
114             }
115 6         38 ($ws,\@wheel);
116             }
117              
118             sub factors_wheel
119             {
120 24     24 1 5922 my $n = shift;
121 24   50     68 my $o = abs(shift || 1);
122              
123 24 50       91 $n = Math::BigInt->new($n) unless ref $n;
124 24         1560 my $two = Math::BigInt->new(2);
125 24         1033 my $three = Math::BigInt->new(3);
126              
127 24         1032 my @factors = ();
128 24         53 my $x = $n->copy();
129              
130 24 100       488 return ($x) if $x < 4;
131 8         809 my ($i,$y,$w,$div,$rem);
132              
133             #print "Using a wheel of order $o, length ";
134 8         15 my $wheel = wheel($o);
135             #print scalar @$wheel,":\n";
136 8         193 my ($ws,$add) = _transform_wheel($wheel); undef $wheel;
  8         59  
137 8         14 my $we = scalar @$add - 1;
138              
139             # reduce to odd number (after that, no odd left-over divisior will ocur)
140 8   66     28 while (($x->is_even) && (!$x->is_zero))
141             {
142 4         93 push @factors, $two->copy();
143             #print "factoring $x (",$x->length(),")\n";
144             #print "2\n";
145 4         70 $x /= $two;
146             }
147             # 8 => 6 => 3, 7, 6 => 3, 5, 4 => 2 => 1, 3, 2 => 1, are all prime
148             # so the first number interesting for us is 9
149             #my $op = 0;
150             OUTER:
151 8         642 while ($x > 8)
152             {
153             #print "factoring $x (",$x->length(),")\n";
154 8         841 $i = $three; $w = 0;
  8         12  
155 8         18 while ($i < $x) # should be sqrt()
156             {
157             # $steps++;
158             # $op = 0, print "$i\r" if $op++ == 1024;
159 8         216 $y = $x->copy();
160 8         135 ($div,$rem) = $y->bdiv($i);
161 8 50       1236 if ($rem == 0)
162             {
163             #print "$i\n";
164 8         1248 push @factors,$i;
165 8         19 $x = $div; next OUTER;
  8         27  
166             }
167             #print "$i + ",$add->[$w]," ($w)\n";
168             #$i += 2; # trial div by odd numbers
169 0         0 $i += $add->[$w];
170             #print "restart $w $ws\n" if $w == $we; # wheel of 2,3,5,7...
171 0 0       0 $w = $ws if $w++ == $we; # wheel of 2,3,5,7...
172             #exit if $i > 100000;
173             }
174 0         0 last;
175             }
176 8 50 33     779 push @factors,$x if $x != 1 || $n == 1;
177 8         922 @factors;
178             }
179              
180             sub _factor
181       0     {
182             # later: factor ( n => $n, algorithmn => 'wheel', order => 3 );
183             }
184              
185             1;
186              
187             __END__