File Coverage

blib/lib/bigint.pl
Criterion Covered Total %
statement 148 152 97.3
branch 89 94 94.6
condition 20 21 95.2
subroutine 17 17 100.0
pod 0 6 0.0
total 274 290 94.4


line stmt bran cond sub pod time code
1             package bigint;
2             #
3             # This library is no longer being maintained, and is included for backward
4             # compatibility with Perl 4 programs which may require it.
5             #
6             # In particular, this should not be used as an example of modern Perl
7             # programming techniques.
8             #
9             # Suggested alternative: Math::BigInt
10             #
11             # arbitrary size integer math package
12             #
13             # by Mark Biggar
14             #
15             # Canonical Big integer value are strings of the form
16             # /^[+-]\d+$/ with leading zeros suppressed
17             # Input values to these routines may be strings of the form
18             # /^\s*[+-]?[\d\s]+$/.
19             # Examples:
20             # '+0' canonical zero value
21             # ' -123 123 123' canonical value '-123123123'
22             # '1 23 456 7890' canonical value '+1234567890'
23             # Output values always in canonical form
24             #
25             # Actual math is done in an internal format consisting of an array
26             # whose first element is the sign (/^[+-]$/) and whose remaining
27             # elements are base 100000 digits with the least significant digit first.
28             # The string 'NaN' is used to represent the result when input arguments
29             # are not numbers, as well as the result of dividing by zero
30             #
31             # routines provided are:
32             #
33             # bneg(BINT) return BINT negation
34             # babs(BINT) return BINT absolute value
35             # bcmp(BINT,BINT) return CODE compare numbers (undef,<0,=0,>0)
36             # badd(BINT,BINT) return BINT addition
37             # bsub(BINT,BINT) return BINT subtraction
38             # bmul(BINT,BINT) return BINT multiplication
39             # bdiv(BINT,BINT) return (BINT,BINT) division (quo,rem) just quo if scalar
40             # bmod(BINT,BINT) return BINT modulus
41             # bgcd(BINT,BINT) return BINT greatest common divisor
42             # bnorm(BINT) return BINT normalization
43             #
44              
45             # overcome a floating point problem on certain osnames (posix-bc, os390)
46             BEGIN {
47 3     3   521 my $x = 100000.0;
48 3 50       5210 my $use_mult = int($x*1e-5)*1e5 == $x ? 1 : 0;
49             }
50              
51             $zero = 0;
52              
53            
54             # normalize string form of number. Strip leading zeros. Strip any
55             # white space and add a sign, if missing.
56             # Strings that are not numbers result the value 'NaN'.
57              
58             sub main::bnorm { #(num_str) return num_str
59 2982     2982   14148 local($_) = @_;
60 2982         4126 s/\s+//g; # strip white space
61 2982 100       12363 if (s/^([+-]?)0*(\d+)$/$1$2/) { # test if number
62 2944 100       5230 substr($_,0,0) = '+' unless $1; # Add missing sign
63 2944         3286 s/^-0/+0/;
64 2944         8229 $_;
65             } else {
66 38         113 'NaN';
67             }
68             }
69              
70             # Convert a number from string format to internal base 100000 format.
71             # Assumes normalized value as input.
72             sub internal { #(num_str) return int_num_array
73 1448     1448 0 1854 local($d) = @_;
74 1448         2213 ($is,$il) = (substr($d,0,1),length($d)-2);
75 1448         1823 substr($d,0,1) = '';
76 1448         5482 ($is, reverse(unpack("a" . ($il%5+1) . ("a5" x ($il/5)), $d)));
77             }
78              
79             # Convert a number from internal base 100000 format to string format.
80             # This routine scribbles all over input array.
81             sub external { #(int_num_array) return num_str
82 862     862 0 1052 $es = shift;
83 862   66     4463 grep($_ > 9999 || ($_ = substr('0000'.$_,-5)), @_); # zero pad
84 862         2361 &::bnorm(join('', $es, reverse(@_))); # reverse concat and normalize
85             }
86              
87             # Negate input value.
88             sub main::bneg { #(num_str) return num_str
89 43     43   2701 local($_) = &::bnorm(@_);
90 43 100       141 vec($_,0,8) ^= ord('+') ^ ord('-') unless $_ eq '+0';
91 43 100       146 s/^./N/ unless /^[-+]/; # works both in ASCII and EBCDIC
92 43         103 $_;
93             }
94              
95             # Returns the absolute value of the input.
96             sub main::babs { #(num_str) return num_str
97 106     106   2769 &abs(&::bnorm(@_));
98             }
99              
100             sub abs { # post-normalized abs for internal use
101 600     600 0 796 local($_) = @_;
102 600         793 s/^-/+/;
103 600         1043 $_;
104             }
105            
106             # Compares 2 values. Returns one of undef, <0, =0, >0. (suitable for sort)
107             sub main::bcmp { #(num_str, num_str) return cond_code
108 207     207   10350 local($x,$y) = (&::bnorm($_[0]),&::bnorm($_[1]));
109 207 100       396 if ($x eq 'NaN') {
    100          
110 2         12 undef;
111             } elsif ($y eq 'NaN') {
112 1         7 undef;
113             } else {
114 204         289 &cmp($x,$y);
115             }
116             }
117              
118             sub cmp { # post-normalized compare for internal use
119 466     466 0 715 local($cx, $cy) = @_;
120 466 100       894 return 0 if ($cx eq $cy);
121              
122 404         781 local($sx, $sy) = (substr($cx, 0, 1), substr($cy, 0, 1));
123 404         419 local($ld);
124              
125 404 100       571 if ($sx eq '+') {
126 394 100 100     1027 return 1 if ($sy eq '-' || $cy eq '+0');
127 381         483 $ld = length($cx) - length($cy);
128 381 100       972 return $ld if ($ld);
129 138         488 return $cx cmp $cy;
130             } else { # $sx eq '-'
131 10 100       28 return -1 if ($sy eq '+');
132 8         16 $ld = length($cy) - length($cx);
133 8 100       64 return $ld if ($ld);
134 4         42 return $cy cmp $cx;
135             }
136              
137             }
138              
139             sub main::badd { #(num_str, num_str) return num_str
140 266     266   17599 local(*x, *y); ($x, $y) = (&::bnorm($_[0]),&::bnorm($_[1]));
  266         419  
141 266 100       523 if ($x eq 'NaN') {
    100          
142 4         29 'NaN';
143             } elsif ($y eq 'NaN') {
144 2         15 'NaN';
145             } else {
146 260         407 @x = &internal($x); # convert to internal form
147 260         426 @y = &internal($y);
148 260         465 local($sx, $sy) = (shift @x, shift @y); # get signs
149 260 100       423 if ($sx eq $sy) {
150 189         298 &external($sx, &add(*x, *y)); # if same sign add
151             } else {
152 71         114 ($x, $y) = (&abs($x),&abs($y)); # make abs
153 71 100       114 if (&cmp($y,$x) > 0) {
154 34         63 &external($sy, &sub(*y, *x));
155             } else {
156 37         67 &external($sx, &sub(*x, *y));
157             }
158             }
159             }
160             }
161              
162             sub main::bsub { #(num_str, num_str) return num_str
163 37     37   17488 &::badd($_[0],&::bneg($_[1]));
164             }
165              
166             # GCD -- Euclids algorithm Knuth Vol 2 pg 296
167             sub main::bgcd { #(num_str, num_str) return num_str
168 11     11   5195 local($x,$y) = (&::bnorm($_[0]),&::bnorm($_[1]));
169 11 100 100     41 if ($x eq 'NaN' || $y eq 'NaN') {
170 3         21 'NaN';
171             } else {
172 8         21 ($x, $y) = ($y,&::bmod($x,$y)) while $y ne '+0';
173 8         60 $x;
174             }
175             }
176            
177             # routine to add two base 1e5 numbers
178             # stolen from Knuth Vol 2 Algorithm A pg 231
179             # there are separate routines to add and sub as per Kunth pg 233
180             sub add { #(int_num_array, int_num_array) return int_num_array
181 189     189 0 272 local(*x, *y) = @_;
182 189         296 $car = 0;
183 189         292 for $x (@x) {
184 501 100 100     792 last unless @y || $car;
185 453 100       1113 $x -= 1e5 if $car = (($x += shift(@y) + $car) >= 1e5) ? 1 : 0;
    100          
186             }
187 189         306 for $y (@y) {
188 13 100       24 last unless $car;
189 7 100       17 $y -= 1e5 if $car = (($y += $car) >= 1e5) ? 1 : 0;
    100          
190             }
191 189         477 (@x, @y, $car);
192             }
193              
194             # subtract base 1e5 numbers -- stolen from Knuth Vol 2 pg 232, $x > $y
195             sub sub { #(int_num_array, int_num_array) return int_num_array
196 71     71 0 141 local(*sx, *sy) = @_;
197 71         89 $bar = 0;
198 71         111 for $sx (@sx) {
199 106 100 100     190 last unless @y || $bar;
200 100 100       302 $sx += 1e5 if $bar = (($sx -= shift(@sy) + $bar) < 0);
201             }
202 71         156 @sx;
203             }
204              
205             # multiply two numbers -- stolen from Knuth Vol 2 pg 233
206             sub main::bmul { #(num_str, num_str) return num_str
207 302     302   15433 local(*x, *y); ($x, $y) = (&::bnorm($_[0]), &::bnorm($_[1]));
  302         477  
208 302 100       566 if ($x eq 'NaN') {
    100          
209 2         16 'NaN';
210             } elsif ($y eq 'NaN') {
211 1         8 'NaN';
212             } else {
213 299         439 @x = &internal($x);
214 299         453 @y = &internal($y);
215 299 100       603 local($signr) = (shift @x ne shift @y) ? '-' : '+';
216 299         379 @prod = ();
217 299         488 for $x (@x) {
218 715         790 ($car, $cty) = (0, 0);
219 715         794 for $y (@y) {
220 741         1260 $prod = $x * $y + $prod[$cty] + $car;
221 741 50       923 if ($use_mult) {
222 0         0 $prod[$cty++] =
223             $prod - ($car = int($prod * 1e-5)) * 1e5;
224             }
225             else {
226 741         1320 $prod[$cty++] =
227             $prod - ($car = int($prod / 1e5)) * 1e5;
228             }
229             }
230 715 100       1051 $prod[$cty] += $car if $car;
231 715         871 $x = shift @prod;
232             }
233 299         481 &external($signr, @x, @prod);
234             }
235             }
236              
237             # modulus
238             sub main::bmod { #(num_str, num_str) return num_str
239 52     52   16781 (&::bdiv(@_))[1];
240             }
241            
242             sub main::bdiv { #(dividend: num_str, divisor: num_str) return num_str
243 188     188   17236 local (*x, *y); ($x, $y) = (&::bnorm($_[0]), &::bnorm($_[1]));
  188         335  
244 188 100 100     821 return wantarray ? ('NaN','NaN') : 'NaN'
    100 100        
245             if ($x eq 'NaN' || $y eq 'NaN' || $y eq '+0');
246 176 100       284 return wantarray ? ('+0',$x) : '+0' if (&cmp(&abs($x),&abs($y)) < 0);
    100          
247 165         283 @x = &internal($x); @y = &internal($y);
  165         277  
248 165         229 $srem = $y[0];
249 165 100       303 $sr = (shift @x ne shift @y) ? '-' : '+';
250 165         272 $car = $bar = $prd = 0;
251 165 100       462 if (($dd = int(1e5/($y[$#y]+1))) != 1) {
252 162         268 for $x (@x) {
253 1008         1159 $x = $x * $dd + $car;
254 1008 50       1125 if ($use_mult) {
255 0         0 $x -= ($car = int($x * 1e-5)) * 1e5;
256             }
257             else {
258 1008         1390 $x -= ($car = int($x / 1e5)) * 1e5;
259             }
260             }
261 162         281 push(@x, $car); $car = 0;
  162         164  
262 162         208 for $y (@y) {
263 288         316 $y = $y * $dd + $car;
264 288 50       324 if ($use_mult) {
265 0         0 $y -= ($car = int($y * 1e-5)) * 1e5;
266             }
267             else {
268 288         428 $y -= ($car = int($y / 1e5)) * 1e5;
269             }
270             }
271             }
272             else {
273 3         9 push(@x, 0);
274             }
275 165         293 @q = (); ($v2,$v1) = @y[-2,-1];
  165         244  
276 165         278 while ($#x > $#y) {
277 897         1252 ($u2,$u1,$u0) = @x[-3..-1];
278 897 100       1427 $q = (($u0 == $v1) ? 99999 : int(($u0*1e5+$u1)/$v1));
279 897         1600 --$q while ($v2*$q > ($u0*1e5+$u1-$q*$v1)*1e5+$u2);
280 897 100       1127 if ($q) {
281 570         618 ($car, $bar) = (0,0);
282 570         1029 for ($y = 0, $x = $#x-$#y-1; $y <= $#y; ++$y,++$x) {
283 2583         2756 $prd = $q * $y[$y] + $car;
284 2583 50       2817 if ($use_mult) {
285 0         0 $prd -= ($car = int($prd * 1e-5)) * 1e5;
286             }
287             else {
288 2583         2952 $prd -= ($car = int($prd / 1e5)) * 1e5;
289             }
290 2583 100       5345 $x[$x] += 1e5 if ($bar = (($x[$x] -= $prd + $bar) < 0));
291             }
292 570 100       877 if ($x[$#x] < $car + $bar) {
293 5         7 $car = 0; --$q;
  5         7  
294 5         15 for ($y = 0, $x = $#x-$#y-1; $y <= $#y; ++$y,++$x) {
295 44 100       84 $x[$x] -= 1e5
296             if ($car = (($x[$x] += $y[$y] + $car) > 1e5));
297             }
298             }
299             }
300 897         936 pop(@x); unshift(@q, $q);
  897         1649  
301             }
302 165 100       227 if (wantarray) {
303 138         212 @d = ();
304 138 100       174 if ($dd != 1) {
305 136         148 $car = 0;
306 136         187 for $x (reverse @x) {
307 262         293 $prd = $car * 1e5 + $x;
308 262         316 $car = $prd - ($tmp = int($prd / $dd)) * $dd;
309 262         345 unshift(@d, $tmp);
310             }
311             }
312             else {
313 2         6 @d = @x;
314             }
315 138         245 (&external($sr, @q), &external($srem, @d, $zero));
316             } else {
317 27         42 &external($sr, @q);
318             }
319             }
320             1;