[win32] merge change#896 from maintbranch
[p5sagit/p5-mst-13.2.git] / lib / Math / BigInt.pm
CommitLineData
a0d0e21e 1package Math::BigInt;
2
a5f75d66 3use overload
748a9306 4'+' => sub {new Math::BigInt &badd},
5'-' => sub {new Math::BigInt
a0d0e21e 6 $_[2]? bsub($_[1],${$_[0]}) : bsub(${$_[0]},$_[1])},
748a9306 7'<=>' => sub {new Math::BigInt
a0d0e21e 8 $_[2]? bcmp($_[1],${$_[0]}) : bcmp(${$_[0]},$_[1])},
748a9306 9'cmp' => sub {new Math::BigInt
a0d0e21e 10 $_[2]? ($_[1] cmp ${$_[0]}) : (${$_[0]} cmp $_[1])},
748a9306 11'*' => sub {new Math::BigInt &bmul},
12'/' => sub {new Math::BigInt
a0d0e21e 13 $_[2]? scalar bdiv($_[1],${$_[0]}) :
14 scalar bdiv(${$_[0]},$_[1])},
748a9306 15'%' => sub {new Math::BigInt
a0d0e21e 16 $_[2]? bmod($_[1],${$_[0]}) : bmod(${$_[0]},$_[1])},
748a9306 17'**' => sub {new Math::BigInt
a0d0e21e 18 $_[2]? bpow($_[1],${$_[0]}) : bpow(${$_[0]},$_[1])},
748a9306 19'neg' => sub {new Math::BigInt &bneg},
20'abs' => sub {new Math::BigInt &babs},
a0d0e21e 21
22qw(
23"" stringify
240+ numify) # Order of arguments unsignificant
a5f75d66 25;
a0d0e21e 26
748a9306 27$NaNOK=1;
28
a0d0e21e 29sub new {
a5f75d66 30 my($class) = shift;
31 my($foo) = bnorm(shift);
748a9306 32 die "Not a number initialized to Math::BigInt" if !$NaNOK && $foo eq "NaN";
a5f75d66 33 bless \$foo, $class;
a0d0e21e 34}
35sub stringify { "${$_[0]}" }
36sub numify { 0 + "${$_[0]}" } # Not needed, additional overhead
37 # comparing to direct compilation based on
38 # stringify
39
a0d0e21e 40$zero = 0;
41
a5f75d66 42
a0d0e21e 43# normalize string form of number. Strip leading zeros. Strip any
44# white space and add a sign, if missing.
45# Strings that are not numbers result the value 'NaN'.
46
47sub bnorm { #(num_str) return num_str
48 local($_) = @_;
49 s/\s+//g; # strip white space
50 if (s/^([+-]?)0*(\d+)$/$1$2/) { # test if number
51 substr($_,$[,0) = '+' unless $1; # Add missing sign
52 s/^-0/+0/;
53 $_;
54 } else {
55 'NaN';
56 }
57}
58
59# Convert a number from string format to internal base 100000 format.
60# Assumes normalized value as input.
61sub internal { #(num_str) return int_num_array
62 local($d) = @_;
63 ($is,$il) = (substr($d,$[,1),length($d)-2);
64 substr($d,$[,1) = '';
65 ($is, reverse(unpack("a" . ($il%5+1) . ("a5" x ($il/5)), $d)));
66}
67
68# Convert a number from internal base 100000 format to string format.
69# This routine scribbles all over input array.
70sub external { #(int_num_array) return num_str
71 $es = shift;
72 grep($_ > 9999 || ($_ = substr('0000'.$_,-5)), @_); # zero pad
73 &bnorm(join('', $es, reverse(@_))); # reverse concat and normalize
74}
75
76# Negate input value.
77sub bneg { #(num_str) return num_str
78 local($_) = &bnorm(@_);
79 vec($_,0,8) ^= ord('+') ^ ord('-') unless $_ eq '+0';
80 s/^H/N/;
81 $_;
82}
83
84# Returns the absolute value of the input.
85sub babs { #(num_str) return num_str
86 &abs(&bnorm(@_));
87}
88
89sub abs { # post-normalized abs for internal use
90 local($_) = @_;
91 s/^-/+/;
92 $_;
93}
a5f75d66 94
a0d0e21e 95# Compares 2 values. Returns one of undef, <0, =0, >0. (suitable for sort)
96sub bcmp { #(num_str, num_str) return cond_code
97 local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
98 if ($x eq 'NaN') {
99 undef;
100 } elsif ($y eq 'NaN') {
101 undef;
102 } else {
103 &cmp($x,$y);
104 }
105}
106
107sub cmp { # post-normalized compare for internal use
108 local($cx, $cy) = @_;
1e2e1ae8 109
110 return 0 if ($cx eq $cy);
111
112 local($sx, $sy) = (substr($cx, 0, 1), substr($cy, 0, 1));
113 local($ld);
114
115 if ($sx eq '+') {
116 return 1 if ($sy eq '-' || $cy eq '+0');
117 $ld = length($cx) - length($cy);
118 return $ld if ($ld);
119 return $cx cmp $cy;
120 } else { # $sx eq '-'
121 return -1 if ($sy eq '+');
122 $ld = length($cy) - length($cx);
123 return $ld if ($ld);
124 return $cy cmp $cx;
125 }
a0d0e21e 126}
127
128sub badd { #(num_str, num_str) return num_str
129 local(*x, *y); ($x, $y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
130 if ($x eq 'NaN') {
131 'NaN';
132 } elsif ($y eq 'NaN') {
133 'NaN';
134 } else {
135 @x = &internal($x); # convert to internal form
136 @y = &internal($y);
137 local($sx, $sy) = (shift @x, shift @y); # get signs
138 if ($sx eq $sy) {
139 &external($sx, &add(*x, *y)); # if same sign add
140 } else {
141 ($x, $y) = (&abs($x),&abs($y)); # make abs
142 if (&cmp($y,$x) > 0) {
143 &external($sy, &sub(*y, *x));
144 } else {
145 &external($sx, &sub(*x, *y));
146 }
147 }
148 }
149}
150
151sub bsub { #(num_str, num_str) return num_str
152 &badd($_[$[],&bneg($_[$[+1]));
153}
154
155# GCD -- Euclids algorithm Knuth Vol 2 pg 296
156sub bgcd { #(num_str, num_str) return num_str
157 local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
158 if ($x eq 'NaN' || $y eq 'NaN') {
159 'NaN';
160 } else {
161 ($x, $y) = ($y,&bmod($x,$y)) while $y ne '+0';
162 $x;
163 }
164}
a5f75d66 165
a0d0e21e 166# routine to add two base 1e5 numbers
167# stolen from Knuth Vol 2 Algorithm A pg 231
168# there are separate routines to add and sub as per Kunth pg 233
169sub add { #(int_num_array, int_num_array) return int_num_array
170 local(*x, *y) = @_;
171 $car = 0;
172 for $x (@x) {
173 last unless @y || $car;
20408e3c 174 $x -= 1e5 if $car = (($x += (@y ? shift(@y) : 0) + $car) >= 1e5) ? 1 : 0;
a0d0e21e 175 }
176 for $y (@y) {
177 last unless $car;
55497cff 178 $y -= 1e5 if $car = (($y += $car) >= 1e5) ? 1 : 0;
a0d0e21e 179 }
180 (@x, @y, $car);
181}
182
183# subtract base 1e5 numbers -- stolen from Knuth Vol 2 pg 232, $x > $y
184sub sub { #(int_num_array, int_num_array) return int_num_array
185 local(*sx, *sy) = @_;
186 $bar = 0;
187 for $sx (@sx) {
20408e3c 188 last unless @sy || $bar;
189 $sx += 1e5 if $bar = (($sx -= (@sy ? shift(@sy) : 0) + $bar) < 0);
a0d0e21e 190 }
191 @sx;
192}
193
194# multiply two numbers -- stolen from Knuth Vol 2 pg 233
195sub bmul { #(num_str, num_str) return num_str
196 local(*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
197 if ($x eq 'NaN') {
198 'NaN';
199 } elsif ($y eq 'NaN') {
200 'NaN';
201 } else {
202 @x = &internal($x);
203 @y = &internal($y);
204 &external(&mul(*x,*y));
205 }
206}
207
208# multiply two numbers in internal representation
209# destroys the arguments, supposes that two arguments are different
210sub mul { #(*int_num_array, *int_num_array) return int_num_array
211 local(*x, *y) = (shift, shift);
212 local($signr) = (shift @x ne shift @y) ? '-' : '+';
213 @prod = ();
214 for $x (@x) {
215 ($car, $cty) = (0, $[);
216 for $y (@y) {
8b6a6e55 217 $prod = $x * $y + ($prod[$cty] || 0) + $car;
a0d0e21e 218 $prod[$cty++] =
219 $prod - ($car = int($prod * 1e-5)) * 1e5;
220 }
221 $prod[$cty] += $car if $car;
222 $x = shift @prod;
223 }
224 ($signr, @x, @prod);
225}
226
227# modulus
228sub bmod { #(num_str, num_str) return num_str
229 (&bdiv(@_))[$[+1];
230}
a5f75d66 231
a0d0e21e 232sub bdiv { #(dividend: num_str, divisor: num_str) return num_str
233 local (*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
234 return wantarray ? ('NaN','NaN') : 'NaN'
235 if ($x eq 'NaN' || $y eq 'NaN' || $y eq '+0');
236 return wantarray ? ('+0',$x) : '+0' if (&cmp(&abs($x),&abs($y)) < 0);
237 @x = &internal($x); @y = &internal($y);
238 $srem = $y[$[];
239 $sr = (shift @x ne shift @y) ? '-' : '+';
240 $car = $bar = $prd = 0;
241 if (($dd = int(1e5/($y[$#y]+1))) != 1) {
242 for $x (@x) {
243 $x = $x * $dd + $car;
244 $x -= ($car = int($x * 1e-5)) * 1e5;
245 }
246 push(@x, $car); $car = 0;
247 for $y (@y) {
248 $y = $y * $dd + $car;
249 $y -= ($car = int($y * 1e-5)) * 1e5;
250 }
251 }
252 else {
253 push(@x, 0);
254 }
255 @q = (); ($v2,$v1) = @y[-2,-1];
256 while ($#x > $#y) {
257 ($u2,$u1,$u0) = @x[-3..-1];
258 $q = (($u0 == $v1) ? 99999 : int(($u0*1e5+$u1)/$v1));
259 --$q while ($v2*$q > ($u0*1e5+$u1-$q*$v1)*1e5+$u2);
260 if ($q) {
261 ($car, $bar) = (0,0);
262 for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) {
263 $prd = $q * $y[$y] + $car;
264 $prd -= ($car = int($prd * 1e-5)) * 1e5;
265 $x[$x] += 1e5 if ($bar = (($x[$x] -= $prd + $bar) < 0));
266 }
267 if ($x[$#x] < $car + $bar) {
268 $car = 0; --$q;
269 for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) {
270 $x[$x] -= 1e5
271 if ($car = (($x[$x] += $y[$y] + $car) > 1e5));
272 }
273 }
274 }
275 pop(@x); unshift(@q, $q);
276 }
277 if (wantarray) {
278 @d = ();
279 if ($dd != 1) {
280 $car = 0;
281 for $x (reverse @x) {
282 $prd = $car * 1e5 + $x;
283 $car = $prd - ($tmp = int($prd / $dd)) * $dd;
284 unshift(@d, $tmp);
285 }
286 }
287 else {
288 @d = @x;
289 }
290 (&external($sr, @q), &external($srem, @d, $zero));
291 } else {
292 &external($sr, @q);
293 }
294}
295
296# compute power of two numbers -- stolen from Knuth Vol 2 pg 233
297sub bpow { #(num_str, num_str) return num_str
298 local(*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
299 if ($x eq 'NaN') {
300 'NaN';
301 } elsif ($y eq 'NaN') {
302 'NaN';
303 } elsif ($x eq '+1') {
304 '+1';
305 } elsif ($x eq '-1') {
306 &bmod($x,2) ? '-1': '+1';
307 } elsif ($y =~ /^-/) {
308 'NaN';
309 } elsif ($x eq '+0' && $y eq '+0') {
310 'NaN';
311 } else {
312 @x = &internal($x);
313 local(@pow2)=@x;
314 local(@pow)=&internal("+1");
315 local($y1,$res,@tmp1,@tmp2)=(1); # need tmp to send to mul
316 while ($y ne '+0') {
317 ($y,$res)=&bdiv($y,2);
318 if ($res ne '+0') {@tmp=@pow2; @pow=&mul(*pow,*tmp);}
319 if ($y ne '+0') {@tmp=@pow2;@pow2=&mul(*pow2,*tmp);}
320 }
321 &external(@pow);
322 }
323}
324
3251;
a5f75d66 326__END__
327
328=head1 NAME
329
330Math::BigInt - Arbitrary size integer math package
331
332=head1 SYNOPSIS
333
334 use Math::BigInt;
335 $i = Math::BigInt->new($string);
336
337 $i->bneg return BINT negation
338 $i->babs return BINT absolute value
339 $i->bcmp(BINT) return CODE compare numbers (undef,<0,=0,>0)
340 $i->badd(BINT) return BINT addition
341 $i->bsub(BINT) return BINT subtraction
342 $i->bmul(BINT) return BINT multiplication
343 $i->bdiv(BINT) return (BINT,BINT) division (quo,rem) just quo if scalar
344 $i->bmod(BINT) return BINT modulus
345 $i->bgcd(BINT) return BINT greatest common divisor
346 $i->bnorm return BINT normalization
347
348=head1 DESCRIPTION
349
350All basic math operations are overloaded if you declare your big
351integers as
352
353 $i = new Math::BigInt '123 456 789 123 456 789';
354
355
356=over 2
357
358=item Canonical notation
359
360Big integer value are strings of the form C</^[+-]\d+$/> with leading
361zeros suppressed.
362
363=item Input
364
365Input values to these routines may be strings of the form
366C</^\s*[+-]?[\d\s]+$/>.
367
368=item Output
369
370Output values always always in canonical form
371
372=back
373
374Actual math is done in an internal format consisting of an array
375whose first element is the sign (/^[+-]$/) and whose remaining
376elements are base 100000 digits with the least significant digit first.
377The string 'NaN' is used to represent the result when input arguments
378are not numbers, as well as the result of dividing by zero.
379
380=head1 EXAMPLES
381
382 '+0' canonical zero value
383 ' -123 123 123' canonical value '-123123123'
384 '1 23 456 7890' canonical value '+1234567890'
385
386
387=head1 BUGS
388
389The current version of this module is a preliminary version of the
390real thing that is currently (as of perl5.002) under development.
391
392=head1 AUTHOR
393
394Mark Biggar, overloaded interface by Ilya Zakharevich.
395
396=cut