[p5sagit/p5-mst-13.2.git] / lib / Math / BigInt.pm

package Math::BigInt;
$VERSION='0.01';

use overload
'+'	=>	sub {new Math::BigInt &badd},
'-'	=>	sub {new Math::BigInt
		       $_[2]? bsub($_[1],${$_[0]}) : bsub(${$_[0]},$_[1])},
'<=>'	=>	sub {$_[2]? bcmp($_[1],${$_[0]}) : bcmp(${$_[0]},$_[1])},
'cmp'	=>	sub {$_[2]? ($_[1] cmp ${$_[0]}) : (${$_[0]} cmp $_[1])},
'*'	=>	sub {new Math::BigInt &bmul},
'/'	=>	sub {new Math::BigInt 
		       $_[2]? scalar bdiv($_[1],${$_[0]}) :
			 scalar bdiv(${$_[0]},$_[1])},
'%'	=>	sub {new Math::BigInt
		       $_[2]? bmod($_[1],${$_[0]}) : bmod(${$_[0]},$_[1])},
'**'	=>	sub {new Math::BigInt
		       $_[2]? bpow($_[1],${$_[0]}) : bpow(${$_[0]},$_[1])},
'neg'	=>	sub {new Math::BigInt &bneg},
'abs'	=>	sub {new Math::BigInt &babs},
'<<'	=>	sub {new Math::BigInt
		       $_[2]? blsft($_[1],${$_[0]}) : blsft(${$_[0]},$_[1])},
'>>'	=>	sub {new Math::BigInt
		       $_[2]? brsft($_[1],${$_[0]}) : brsft(${$_[0]},$_[1])},
'&'	=>	sub {new Math::BigInt &band},
'|'	=>	sub {new Math::BigInt &bior},
'^'	=>	sub {new Math::BigInt &bxor},
'~'	=>	sub {new Math::BigInt &bnot},
'int'	=>	sub { shift },

qw(
""	stringify
0+	numify)			# Order of arguments unsignificant
;

$NaNOK=1;

sub new {
  my($class) = shift;
  my($foo) = bnorm(shift);
  die "Not a number initialized to Math::BigInt" if !$NaNOK && $foo eq "NaN";
  bless \$foo, $class;
}
sub stringify { "${$_[0]}" }
sub numify { 0 + "${$_[0]}" }	# Not needed, additional overhead
				# comparing to direct compilation based on
				# stringify
sub import {
  shift;
  return unless @_;
  die "unknown import: @_" unless @_ == 1 and $_[0] eq ':constant';
  overload::constant integer => sub {Math::BigInt->new(shift)};
}

$zero = 0;

# overcome a floating point problem on certain osnames (posix-bc, os390)
BEGIN {
    my $x = 100000.0;
    my $use_mult = int($x*1e-5)*1e5 == $x ? 1 : 0;
}

# normalize string form of number.   Strip leading zeros.  Strip any
#   white space and add a sign, if missing.
# Strings that are not numbers result the value 'NaN'.

sub bnorm { #(num_str) return num_str
    local($_) = @_;
    s/\s+//g;                           # strip white space
    if (s/^([+-]?)0*(\d+)$/$1$2/) {     # test if number
	substr($_,$[,0) = '+' unless $1; # Add missing sign
	s/^-0/+0/;
	$_;
    } else {
	'NaN';
    }
}

# Convert a number from string format to internal base 100000 format.
#   Assumes normalized value as input.
sub internal { #(num_str) return int_num_array
    local($d) = @_;
    ($is,$il) = (substr($d,$[,1),length($d)-2);
    substr($d,$[,1) = '';
    ($is, reverse(unpack("a" . ($il%5+1) . ("a5" x ($il/5)), $d)));
}

# Convert a number from internal base 100000 format to string format.
#   This routine scribbles all over input array.
sub external { #(int_num_array) return num_str
    $es = shift;
    grep($_ > 9999 || ($_ = substr('0000'.$_,-5)), @_);   # zero pad
    &bnorm(join('', $es, reverse(@_)));    # reverse concat and normalize
}

# Negate input value.
sub bneg { #(num_str) return num_str
    local($_) = &bnorm(@_);
    return $_ if $_ eq '+0' or $_ eq 'NaN';
    vec($_,0,8) ^= ord('+') ^ ord('-');
    $_;
}

# Returns the absolute value of the input.
sub babs { #(num_str) return num_str
    &abs(&bnorm(@_));
}

sub abs { # post-normalized abs for internal use
    local($_) = @_;
    s/^-/+/;
    $_;
}

# Compares 2 values.  Returns one of undef, <0, =0, >0. (suitable for sort)
sub bcmp { #(num_str, num_str) return cond_code
    local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
    if ($x eq 'NaN') {
	undef;
    } elsif ($y eq 'NaN') {
	undef;
    } else {
	&cmp($x,$y) <=> 0;
    }
}

sub cmp { # post-normalized compare for internal use
    local($cx, $cy) = @_;
    
    return 0 if ($cx eq $cy);

    local($sx, $sy) = (substr($cx, 0, 1), substr($cy, 0, 1));
    local($ld);

    if ($sx eq '+') {
      return  1 if ($sy eq '-' || $cy eq '+0');
      $ld = length($cx) - length($cy);
      return $ld if ($ld);
      return $cx cmp $cy;
    } else { # $sx eq '-'
      return -1 if ($sy eq '+');
      $ld = length($cy) - length($cx);
      return $ld if ($ld);
      return $cy cmp $cx;
    }
}

sub badd { #(num_str, num_str) return num_str
    local(*x, *y); ($x, $y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
    if ($x eq 'NaN') {
	'NaN';
    } elsif ($y eq 'NaN') {
	'NaN';
    } else {
	@x = &internal($x);             # convert to internal form
	@y = &internal($y);
	local($sx, $sy) = (shift @x, shift @y); # get signs
	if ($sx eq $sy) {
	    &external($sx, &add(*x, *y)); # if same sign add
	} else {
	    ($x, $y) = (&abs($x),&abs($y)); # make abs
	    if (&cmp($y,$x) > 0) {
		&external($sy, &sub(*y, *x));
	    } else {
		&external($sx, &sub(*x, *y));
	    }
	}
    }
}

sub bsub { #(num_str, num_str) return num_str
    &badd($_[$[],&bneg($_[$[+1]));    
}

# GCD -- Euclids algorithm Knuth Vol 2 pg 296
sub bgcd { #(num_str, num_str) return num_str
    local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
    if ($x eq 'NaN' || $y eq 'NaN') {
	'NaN';
    } else {
	($x, $y) = ($y,&bmod($x,$y)) while $y ne '+0';
	$x;
    }
}

# routine to add two base 1e5 numbers
#   stolen from Knuth Vol 2 Algorithm A pg 231
#   there are separate routines to add and sub as per Kunth pg 233
sub add { #(int_num_array, int_num_array) return int_num_array
    local(*x, *y) = @_;
    $car = 0;
    for $x (@x) {
	last unless @y || $car;
	$x -= 1e5 if $car = (($x += (@y ? shift(@y) : 0) + $car) >= 1e5) ? 1 : 0;
    }
    for $y (@y) {
	last unless $car;
	$y -= 1e5 if $car = (($y += $car) >= 1e5) ? 1 : 0;
    }
    (@x, @y, $car);
}

# subtract base 1e5 numbers -- stolen from Knuth Vol 2 pg 232, $x > $y
sub sub { #(int_num_array, int_num_array) return int_num_array
    local(*sx, *sy) = @_;
    $bar = 0;
    for $sx (@sx) {
	last unless @sy || $bar;
	$sx += 1e5 if $bar = (($sx -= (@sy ? shift(@sy) : 0) + $bar) < 0);
    }
    @sx;
}

# multiply two numbers -- stolen from Knuth Vol 2 pg 233
sub bmul { #(num_str, num_str) return num_str
    local(*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
    if ($x eq 'NaN') {
	'NaN';
    } elsif ($y eq 'NaN') {
	'NaN';
    } else {
	@x = &internal($x);
	@y = &internal($y);
	&external(&mul(*x,*y));
    }
}

# multiply two numbers in internal representation
# destroys the arguments, supposes that two arguments are different
sub mul { #(*int_num_array, *int_num_array) return int_num_array
    local(*x, *y) = (shift, shift);
    local($signr) = (shift @x ne shift @y) ? '-' : '+';
    @prod = ();
    for $x (@x) {
      ($car, $cty) = (0, $[);
      for $y (@y) {
	$prod = $x * $y + ($prod[$cty] || 0) + $car;
        if ($use_mult) {
	$prod[$cty++] =
	  $prod - ($car = int($prod * 1e-5)) * 1e5;
        }
        else {
	$prod[$cty++] =
	  $prod - ($car = int($prod / 1e5)) * 1e5;
        }
      }
      $prod[$cty] += $car if $car;
      $x = shift @prod;
    }
    ($signr, @x, @prod);
}

# modulus
sub bmod { #(num_str, num_str) return num_str
    (&bdiv(@_))[$[+1];
}

sub bdiv { #(dividend: num_str, divisor: num_str) return num_str
    local (*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
    return wantarray ? ('NaN','NaN') : 'NaN'
	if ($x eq 'NaN' || $y eq 'NaN' || $y eq '+0');
    return wantarray ? ('+0',$x) : '+0' if (&cmp(&abs($x),&abs($y)) < 0);
    @x = &internal($x); @y = &internal($y);
    $srem = $y[$[];
    $sr = (shift @x ne shift @y) ? '-' : '+';
    $car = $bar = $prd = 0;
    if (($dd = int(1e5/($y[$#y]+1))) != 1) {
	for $x (@x) {
	    $x = $x * $dd + $car;
            if ($use_mult) {
	    $x -= ($car = int($x * 1e-5)) * 1e5;
            }
            else {
	    $x -= ($car = int($x / 1e5)) * 1e5;
            }
	}
	push(@x, $car); $car = 0;
	for $y (@y) {
	    $y = $y * $dd + $car;
            if ($use_mult) {
	    $y -= ($car = int($y * 1e-5)) * 1e5;
            }
            else {
	    $y -= ($car = int($y / 1e5)) * 1e5;
            }
	}
    }
    else {
	push(@x, 0);
    }
    @q = (); ($v2,$v1) = @y[-2,-1];
    $v2 = 0 unless $v2;
    while ($#x > $#y) {
	($u2,$u1,$u0) = @x[-3..-1];
	$u2 = 0 unless $u2;
	$q = (($u0 == $v1) ? 99999 : int(($u0*1e5+$u1)/$v1));
	--$q while ($v2*$q > ($u0*1e5+$u1-$q*$v1)*1e5+$u2);
	if ($q) {
	    ($car, $bar) = (0,0);
	    for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) {
		$prd = $q * $y[$y] + $car;
                if ($use_mult) {
		$prd -= ($car = int($prd * 1e-5)) * 1e5;
                }
                else {
		$prd -= ($car = int($prd / 1e5)) * 1e5;
                }
		$x[$x] += 1e5 if ($bar = (($x[$x] -= $prd + $bar) < 0));
	    }
	    if ($x[$#x] < $car + $bar) {
		$car = 0; --$q;
		for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) {
		    $x[$x] -= 1e5
			if ($car = (($x[$x] += $y[$y] + $car) > 1e5));
		}
	    }   
	}
	pop(@x); unshift(@q, $q);
    }
    if (wantarray) {
	@d = ();
	if ($dd != 1) {
	    $car = 0;
	    for $x (reverse @x) {
		$prd = $car * 1e5 + $x;
		$car = $prd - ($tmp = int($prd / $dd)) * $dd;
		unshift(@d, $tmp);
	    }
	}
	else {
	    @d = @x;
	}
	(&external($sr, @q), &external($srem, @d, $zero));
    } else {
	&external($sr, @q);
    }
}

# compute power of two numbers -- stolen from Knuth Vol 2 pg 233
sub bpow { #(num_str, num_str) return num_str
    local(*x, *y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
    if ($x eq 'NaN') {
	'NaN';
    } elsif ($y eq 'NaN') {
	'NaN';
    } elsif ($x eq '+1') {
	'+1';
    } elsif ($x eq '-1') {
	&bmod($x,2) ? '-1': '+1';
    } elsif ($y =~ /^-/) {
	'NaN';
    } elsif ($x eq '+0' && $y eq '+0') {
	'NaN';
    } else {
	@x = &internal($x);
	local(@pow2)=@x;
	local(@pow)=&internal("+1");
	local($y1,$res,@tmp1,@tmp2)=(1); # need tmp to send to mul
	while ($y ne '+0') {
	  ($y,$res)=&bdiv($y,2);
	  if ($res ne '+0') {@tmp=@pow2; @pow=&mul(*pow,*tmp);}
	  if ($y ne '+0') {@tmp=@pow2;@pow2=&mul(*pow2,*tmp);}
	}
	&external(@pow);
    }
}

# compute x << y, y >= 0
sub blsft { #(num_str, num_str) return num_str
    &bmul($_[$[], &bpow(2, $_[$[+1]));
}

# compute x >> y, y >= 0
sub brsft { #(num_str, num_str) return num_str
    &bdiv($_[$[], &bpow(2, $_[$[+1]));
}

# compute x & y
sub band { #(num_str, num_str) return num_str
    local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
    if ($x eq 'NaN' || $y eq 'NaN') {
	'NaN';
    } else {
	while ($x ne '+0' && $y ne '+0') {
	    ($x, $xr) = &bdiv($x, 0x10000);
	    ($y, $yr) = &bdiv($y, 0x10000);
	    $r = &badd(&bmul(int $xr & $yr, $m), $r);
	    $m = &bmul($m, 0x10000);
	}
	$r;
    }
}

# compute x | y
sub bior { #(num_str, num_str) return num_str
    local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
    if ($x eq 'NaN' || $y eq 'NaN') {
	'NaN';
    } else {
	while ($x ne '+0' || $y ne '+0') {
	    ($x, $xr) = &bdiv($x, 0x10000);
	    ($y, $yr) = &bdiv($y, 0x10000);
	    $r = &badd(&bmul(int $xr | $yr, $m), $r);
	    $m = &bmul($m, 0x10000);
	}
	$r;
    }
}

# compute x ^ y
sub bxor { #(num_str, num_str) return num_str
    local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
    if ($x eq 'NaN' || $y eq 'NaN') {
	'NaN';
    } else {
	while ($x ne '+0' || $y ne '+0') {
	    ($x, $xr) = &bdiv($x, 0x10000);
	    ($y, $yr) = &bdiv($y, 0x10000);
	    $r = &badd(&bmul(int $xr ^ $yr, $m), $r);
	    $m = &bmul($m, 0x10000);
	}
	$r;
    }
}

# represent ~x as twos-complement number
sub bnot { #(num_str) return num_str
    &bsub(-1,$_[$[]);
}

1;
__END__

=head1 NAME

Math::BigInt - Arbitrary size integer math package

=head1 SYNOPSIS

  use Math::BigInt;
  $i = Math::BigInt->new($string);

  $i->bneg return BINT               negation
  $i->babs return BINT               absolute value
  $i->bcmp(BINT) return CODE         compare numbers (undef,<0,=0,>0)
  $i->badd(BINT) return BINT         addition
  $i->bsub(BINT) return BINT         subtraction
  $i->bmul(BINT) return BINT         multiplication
  $i->bdiv(BINT) return (BINT,BINT)  division (quo,rem) just quo if scalar
  $i->bmod(BINT) return BINT         modulus
  $i->bgcd(BINT) return BINT         greatest common divisor
  $i->bnorm return BINT              normalization
  $i->blsft(BINT) return BINT        left shift
  $i->brsft(BINT) return (BINT,BINT) right shift (quo,rem) just quo if scalar
  $i->band(BINT) return BINT         bit-wise and
  $i->bior(BINT) return BINT         bit-wise inclusive or
  $i->bxor(BINT) return BINT         bit-wise exclusive or
  $i->bnot return BINT               bit-wise not

=head1 DESCRIPTION

All basic math operations are overloaded if you declare your big
integers as

  $i = new Math::BigInt '123 456 789 123 456 789';


=over 2

=item Canonical notation

Big integer value are strings of the form C</^[+-]\d+$/> with leading
zeros suppressed.

=item Input

Input values to these routines may be strings of the form
C</^\s*[+-]?[\d\s]+$/>.

=item Output

Output values always always in canonical form

=back

Actual math is done in an internal format consisting of an array
whose first element is the sign (/^[+-]$/) and whose remaining 
elements are base 100000 digits with the least significant digit first.
The string 'NaN' is used to represent the result when input arguments 
are not numbers, as well as the result of dividing by zero.

=head1 EXAMPLES

   '+0'                            canonical zero value
   '   -123 123 123'               canonical value '-123123123'
   '1 23 456 7890'                 canonical value '+1234567890'


=head1 Autocreating constants

After C<use Math::BigInt ':constant'> all the integer decimal constants
in the given scope are converted to C<Math::BigInt>.  This conversion
happens at compile time.

In particular

  perl -MMath::BigInt=:constant -e 'print 2**100'

print the integer value of C<2**100>.  Note that without conversion of 
constants the expression 2**100 will be calculated as floating point number.

=head1 BUGS

The current version of this module is a preliminary version of the
real thing that is currently (as of perl5.002) under development.

=head1 AUTHOR

Mark Biggar, overloaded interface by Ilya Zakharevich.

=cut
Commit	Line	Data
a0d0e21e	1	package Math::BigInt;
3deb277d	2	$VERSION='0.01';
a0d0e21e	3
a5f75d66	4	use overload
748a9306	5	'+' => sub {new Math::BigInt &badd},
748a9306	6	'-' => sub {new Math::BigInt
a0d0e21e	7	$_[2]? bsub($_[1],${$_[0]}) : bsub(${$_[0]},$_[1])},
5d7098d5	8	'<=>' => sub {$_[2]? bcmp($_[1],${$_[0]}) : bcmp(${$_[0]},$_[1])},
5d7098d5	9	'cmp' => sub {$_[2]? ($_[1] cmp ${$_[0]}) : (${$_[0]} cmp $_[1])},
748a9306	10	'*' => sub {new Math::BigInt &bmul},
748a9306	11	'/' => sub {new Math::BigInt
a0d0e21e	12	$_[2]? scalar bdiv($_[1],${$_[0]}) :
a0d0e21e	13	scalar bdiv(${$_[0]},$_[1])},
748a9306	14	'%' => sub {new Math::BigInt
a0d0e21e	15	$_[2]? bmod($_[1],${$_[0]}) : bmod(${$_[0]},$_[1])},
748a9306	16	'**' => sub {new Math::BigInt
a0d0e21e	17	$_[2]? bpow($_[1],${$_[0]}) : bpow(${$_[0]},$_[1])},
748a9306	18	'neg' => sub {new Math::BigInt &bneg},
748a9306	19	'abs' => sub {new Math::BigInt &babs},
e16b8f49	20	'<<' => sub {new Math::BigInt
	21	$_[2]? blsft($_[1],${$_[0]}) : blsft(${$_[0]},$_[1])},
	22	'>>' => sub {new Math::BigInt
	23	$_[2]? brsft($_[1],${$_[0]}) : brsft(${$_[0]},$_[1])},
	24	'&' => sub {new Math::BigInt &band},
	25	'\|' => sub {new Math::BigInt &bior},
	26	'^' => sub {new Math::BigInt &bxor},
	27	'~' => sub {new Math::BigInt &bnot},
f216259d	28	'int' => sub { shift },
a0d0e21e	29
	30	qw(
	31	"" stringify
	32	0+ numify) # Order of arguments unsignificant
a5f75d66	33	;
a0d0e21e	34
748a9306	35	$NaNOK=1;
748a9306	36
a0d0e21e	37	sub new {
a5f75d66	38	my($class) = shift;
a5f75d66	39	my($foo) = bnorm(shift);
748a9306	40	die "Not a number initialized to Math::BigInt" if !$NaNOK && $foo eq "NaN";
a5f75d66	41	bless \$foo, $class;
a0d0e21e	42	}
	43	sub stringify { "${$_[0]}" }
	44	sub numify { 0 + "${$_[0]}" } # Not needed, additional overhead
	45	# comparing to direct compilation based on
	46	# stringify
b3ac6de7	47	sub import {
	48	shift;
	49	return unless @_;
	50	die "unknown import: @_" unless @_ == 1 and $_[0] eq ':constant';
	51	overload::constant integer => sub {Math::BigInt->new(shift)};
	52	}
a0d0e21e	53
a0d0e21e	54	$zero = 0;
a0d0e21e	55
1f45ae4a	56	# overcome a floating point problem on certain osnames (posix-bc, os390)
	57	BEGIN {
	58	my $x = 100000.0;
	59	my $use_mult = int($x1e-5)1e5 == $x ? 1 : 0;
	60	}
a5f75d66	61
a0d0e21e	62	# normalize string form of number. Strip leading zeros. Strip any
	63	# white space and add a sign, if missing.
	64	# Strings that are not numbers result the value 'NaN'.
	65
	66	sub bnorm { #(num_str) return num_str
	67	local($_) = @_;
	68	s/\s+//g; # strip white space
	69	if (s/^([+-]?)0*(\d+)$/$1$2/) { # test if number
	70	substr($_,$[,0) = '+' unless $1; # Add missing sign
	71	s/^-0/+0/;
	72	$_;
	73	} else {
	74	'NaN';
	75	}
	76	}
	77
	78	# Convert a number from string format to internal base 100000 format.
	79	# Assumes normalized value as input.
	80	sub internal { #(num_str) return int_num_array
	81	local($d) = @_;
	82	($is,$il) = (substr($d,$[,1),length($d)-2);
	83	substr($d,$[,1) = '';
	84	($is, reverse(unpack("a" . ($il%5+1) . ("a5" x ($il/5)), $d)));
	85	}
	86
	87	# Convert a number from internal base 100000 format to string format.
	88	# This routine scribbles all over input array.
	89	sub external { #(int_num_array) return num_str
	90	$es = shift;
	91	grep($_ > 9999 \|\| ($_ = substr('0000'.$_,-5)), @_); # zero pad
	92	&bnorm(join('', $es, reverse(@_))); # reverse concat and normalize
	93	}
	94
	95	# Negate input value.
	96	sub bneg { #(num_str) return num_str
	97	local($_) = &bnorm(@_);
e3c7ef20	98	return $_ if $_ eq '+0' or $_ eq 'NaN';
e3c7ef20	99	vec($_,0,8) ^= ord('+') ^ ord('-');
a0d0e21e	100	$_;
	101	}
	102
	103	# Returns the absolute value of the input.
	104	sub babs { #(num_str) return num_str
	105	&abs(&bnorm(@_));
	106	}
	107
	108	sub abs { # post-normalized abs for internal use
	109	local($_) = @_;
	110	s/^-/+/;
	111	$_;
	112	}
a5f75d66	113
a0d0e21e	114	# Compares 2 values. Returns one of undef, <0, =0, >0. (suitable for sort)
	115	sub bcmp { #(num_str, num_str) return cond_code
	116	local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
	117	if ($x eq 'NaN') {
	118	undef;
	119	} elsif ($y eq 'NaN') {
	120	undef;
	121	} else {
e3c7ef20	122	&cmp($x,$y) <=> 0;
a0d0e21e	123	}
	124	}
	125
	126	sub cmp { # post-normalized compare for internal use
	127	local($cx, $cy) = @_;
1e2e1ae8	128
	129	return 0 if ($cx eq $cy);
	130
	131	local($sx, $sy) = (substr($cx, 0, 1), substr($cy, 0, 1));
	132	local($ld);
	133
	134	if ($sx eq '+') {
	135	return 1 if ($sy eq '-' \|\| $cy eq '+0');
	136	$ld = length($cx) - length($cy);
	137	return $ld if ($ld);
	138	return $cx cmp $cy;
	139	} else { # $sx eq '-'
	140	return -1 if ($sy eq '+');
	141	$ld = length($cy) - length($cx);
	142	return $ld if ($ld);
	143	return $cy cmp $cx;
	144	}
a0d0e21e	145	}
	146
	147	sub badd { #(num_str, num_str) return num_str
	148	local(x, y); ($x, $y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
	149	if ($x eq 'NaN') {
	150	'NaN';
	151	} elsif ($y eq 'NaN') {
	152	'NaN';
	153	} else {
	154	@x = &internal($x); # convert to internal form
	155	@y = &internal($y);
	156	local($sx, $sy) = (shift @x, shift @y); # get signs
	157	if ($sx eq $sy) {
	158	&external($sx, &add(x, y)); # if same sign add
	159	} else {
	160	($x, $y) = (&abs($x),&abs($y)); # make abs
	161	if (&cmp($y,$x) > 0) {
	162	&external($sy, &sub(y, x));
	163	} else {
	164	&external($sx, &sub(x, y));
	165	}
	166	}
	167	}
	168	}
	169
	170	sub bsub { #(num_str, num_str) return num_str
	171	&badd($_[$[],&bneg($_[$[+1]));
	172	}
	173
	174	# GCD -- Euclids algorithm Knuth Vol 2 pg 296
	175	sub bgcd { #(num_str, num_str) return num_str
	176	local($x,$y) = (&bnorm($_[$[]),&bnorm($_[$[+1]));
	177	if ($x eq 'NaN' \|\| $y eq 'NaN') {
	178	'NaN';
	179	} else {
	180	($x, $y) = ($y,&bmod($x,$y)) while $y ne '+0';
	181	$x;
	182	}
	183	}
a5f75d66	184
a0d0e21e	185	# routine to add two base 1e5 numbers
	186	# stolen from Knuth Vol 2 Algorithm A pg 231
	187	# there are separate routines to add and sub as per Kunth pg 233
	188	sub add { #(int_num_array, int_num_array) return int_num_array
	189	local(x, y) = @_;
	190	$car = 0;
	191	for $x (@x) {
	192	last unless @y \|\| $car;
20408e3c	193	$x -= 1e5 if $car = (($x += (@y ? shift(@y) : 0) + $car) >= 1e5) ? 1 : 0;
a0d0e21e	194	}
	195	for $y (@y) {
	196	last unless $car;
55497cff	197	$y -= 1e5 if $car = (($y += $car) >= 1e5) ? 1 : 0;
a0d0e21e	198	}
	199	(@x, @y, $car);
	200	}
	201
	202	# subtract base 1e5 numbers -- stolen from Knuth Vol 2 pg 232, $x > $y
	203	sub sub { #(int_num_array, int_num_array) return int_num_array
	204	local(sx, sy) = @_;
	205	$bar = 0;
	206	for $sx (@sx) {
20408e3c	207	last unless @sy \|\| $bar;
20408e3c	208	$sx += 1e5 if $bar = (($sx -= (@sy ? shift(@sy) : 0) + $bar) < 0);
a0d0e21e	209	}
	210	@sx;
	211	}
	212
	213	# multiply two numbers -- stolen from Knuth Vol 2 pg 233
	214	sub bmul { #(num_str, num_str) return num_str
	215	local(x, y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
	216	if ($x eq 'NaN') {
	217	'NaN';
	218	} elsif ($y eq 'NaN') {
	219	'NaN';
	220	} else {
	221	@x = &internal($x);
	222	@y = &internal($y);
	223	&external(&mul(x,y));
	224	}
	225	}
	226
	227	# multiply two numbers in internal representation
	228	# destroys the arguments, supposes that two arguments are different
	229	sub mul { #(int_num_array, int_num_array) return int_num_array
	230	local(x, y) = (shift, shift);
	231	local($signr) = (shift @x ne shift @y) ? '-' : '+';
	232	@prod = ();
	233	for $x (@x) {
	234	($car, $cty) = (0, $[);
	235	for $y (@y) {
8b6a6e55	236	$prod = $x * $y + ($prod[$cty] \|\| 0) + $car;
1f45ae4a	237	if ($use_mult) {
a0d0e21e	238	$prod[$cty++] =
a0d0e21e	239	$prod - ($car = int($prod * 1e-5)) * 1e5;
1f45ae4a	240	}
	241	else {
	242	$prod[$cty++] =
	243	$prod - ($car = int($prod / 1e5)) * 1e5;
	244	}
a0d0e21e	245	}
	246	$prod[$cty] += $car if $car;
	247	$x = shift @prod;
	248	}
	249	($signr, @x, @prod);
	250	}
	251
	252	# modulus
	253	sub bmod { #(num_str, num_str) return num_str
	254	(&bdiv(@_))[$[+1];
	255	}
a5f75d66	256
a0d0e21e	257	sub bdiv { #(dividend: num_str, divisor: num_str) return num_str
	258	local (x, y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
	259	return wantarray ? ('NaN','NaN') : 'NaN'
	260	if ($x eq 'NaN' \|\| $y eq 'NaN' \|\| $y eq '+0');
	261	return wantarray ? ('+0',$x) : '+0' if (&cmp(&abs($x),&abs($y)) < 0);
	262	@x = &internal($x); @y = &internal($y);
	263	$srem = $y[$[];
	264	$sr = (shift @x ne shift @y) ? '-' : '+';
	265	$car = $bar = $prd = 0;
	266	if (($dd = int(1e5/($y[$#y]+1))) != 1) {
	267	for $x (@x) {
	268	$x = $x * $dd + $car;
1f45ae4a	269	if ($use_mult) {
a0d0e21e	270	$x -= ($car = int($x * 1e-5)) * 1e5;
1f45ae4a	271	}
	272	else {
	273	$x -= ($car = int($x / 1e5)) * 1e5;
	274	}
a0d0e21e	275	}
	276	push(@x, $car); $car = 0;
	277	for $y (@y) {
	278	$y = $y * $dd + $car;
1f45ae4a	279	if ($use_mult) {
a0d0e21e	280	$y -= ($car = int($y * 1e-5)) * 1e5;
1f45ae4a	281	}
	282	else {
	283	$y -= ($car = int($y / 1e5)) * 1e5;
	284	}
a0d0e21e	285	}
	286	}
	287	else {
	288	push(@x, 0);
	289	}
5d7098d5	290	@q = (); ($v2,$v1) = @y[-2,-1];
0a6a0d52	291	$v2 = 0 unless $v2;
a0d0e21e	292	while ($#x > $#y) {
5d7098d5	293	($u2,$u1,$u0) = @x[-3..-1];
0a6a0d52	294	$u2 = 0 unless $u2;
a0d0e21e	295	$q = (($u0 == $v1) ? 99999 : int(($u0*1e5+$u1)/$v1));
	296	--$q while ($v2$q > ($u01e5+$u1-$q$v1)1e5+$u2);
	297	if ($q) {
	298	($car, $bar) = (0,0);
	299	for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) {
	300	$prd = $q * $y[$y] + $car;
1f45ae4a	301	if ($use_mult) {
a0d0e21e	302	$prd -= ($car = int($prd * 1e-5)) * 1e5;
1f45ae4a	303	}
	304	else {
	305	$prd -= ($car = int($prd / 1e5)) * 1e5;
	306	}
a0d0e21e	307	$x[$x] += 1e5 if ($bar = (($x[$x] -= $prd + $bar) < 0));
	308	}
	309	if ($x[$#x] < $car + $bar) {
	310	$car = 0; --$q;
	311	for ($y = $[, $x = $#x-$#y+$[-1; $y <= $#y; ++$y,++$x) {
	312	$x[$x] -= 1e5
	313	if ($car = (($x[$x] += $y[$y] + $car) > 1e5));
	314	}
	315	}
	316	}
	317	pop(@x); unshift(@q, $q);
	318	}
	319	if (wantarray) {
	320	@d = ();
	321	if ($dd != 1) {
	322	$car = 0;
	323	for $x (reverse @x) {
	324	$prd = $car * 1e5 + $x;
	325	$car = $prd - ($tmp = int($prd / $dd)) * $dd;
	326	unshift(@d, $tmp);
	327	}
	328	}
	329	else {
	330	@d = @x;
	331	}
	332	(&external($sr, @q), &external($srem, @d, $zero));
	333	} else {
	334	&external($sr, @q);
	335	}
	336	}
	337
	338	# compute power of two numbers -- stolen from Knuth Vol 2 pg 233
	339	sub bpow { #(num_str, num_str) return num_str
	340	local(x, y); ($x, $y) = (&bnorm($_[$[]), &bnorm($_[$[+1]));
	341	if ($x eq 'NaN') {
	342	'NaN';
	343	} elsif ($y eq 'NaN') {
	344	'NaN';
	345	} elsif ($x eq '+1') {
	346	'+1';
	347	} elsif ($x eq '-1') {
	348	&bmod($x,2) ? '-1': '+1';
	349	} elsif ($y =~ /^-/) {
	350	'NaN';
	351	} elsif ($x eq '+0' && $y eq '+0') {
	352	'NaN';
	353	} else {
	354	@x = &internal($x);
	355	local(@pow2)=@x;
	356	local(@pow)=&internal("+1");
	357	local($y1,$res,@tmp1,@tmp2)=(1); # need tmp to send to mul
	358	while ($y ne '+0') {
	359	($y,$res)=&bdiv($y,2);
	360	if ($res ne '+0') {@tmp=@pow2; @pow=&mul(pow,tmp);}
	361	if ($y ne '+0') {@tmp=@pow2;@pow2=&mul(pow2,tmp);}
	362	}
	363	&external(@pow);
	364	}
	365	}
	366
e16b8f49	367	# compute x << y, y >= 0
	368	sub blsft { #(num_str, num_str) return num_str
	369	&bmul($_[$[], &bpow(2, $_[$[+1]));
	370	}
	371
	372	# compute x >> y, y >= 0
	373	sub brsft { #(num_str, num_str) return num_str
	374	&bdiv($_[$[], &bpow(2, $_[$[+1]));
	375	}
	376
	377	# compute x & y
	378	sub band { #(num_str, num_str) return num_str
	379	local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
	380	if ($x eq 'NaN' \|\| $y eq 'NaN') {
	381	'NaN';
	382	} else {
	383	while ($x ne '+0' && $y ne '+0') {
	384	($x, $xr) = &bdiv($x, 0x10000);
	385	($y, $yr) = &bdiv($y, 0x10000);
	386	$r = &badd(&bmul(int $xr & $yr, $m), $r);
	387	$m = &bmul($m, 0x10000);
	388	}
	389	$r;
	390	}
	391	}
	392
	393	# compute x \| y
	394	sub bior { #(num_str, num_str) return num_str
	395	local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
	396	if ($x eq 'NaN' \|\| $y eq 'NaN') {
	397	'NaN';
	398	} else {
	399	while ($x ne '+0' \|\| $y ne '+0') {
	400	($x, $xr) = &bdiv($x, 0x10000);
	401	($y, $yr) = &bdiv($y, 0x10000);
	402	$r = &badd(&bmul(int $xr \| $yr, $m), $r);
	403	$m = &bmul($m, 0x10000);
	404	}
	405	$r;
	406	}
	407	}
	408
	409	# compute x ^ y
	410	sub bxor { #(num_str, num_str) return num_str
	411	local($x,$y,$r,$m,$xr,$yr) = (&bnorm($_[$[]),&bnorm($_[$[+1]),0,1);
	412	if ($x eq 'NaN' \|\| $y eq 'NaN') {
	413	'NaN';
	414	} else {
	415	while ($x ne '+0' \|\| $y ne '+0') {
	416	($x, $xr) = &bdiv($x, 0x10000);
	417	($y, $yr) = &bdiv($y, 0x10000);
	418	$r = &badd(&bmul(int $xr ^ $yr, $m), $r);
	419	$m = &bmul($m, 0x10000);
	420	}
	421	$r;
	422	}
	423	}
	424
	425	# represent ~x as twos-complement number
	426	sub bnot { #(num_str) return num_str
	427	&bsub(-1,$_[$[]);
	428	}
	429
a0d0e21e	430	1;
a5f75d66	431	__END__
	432
	433	=head1 NAME
	434
	435	Math::BigInt - Arbitrary size integer math package
	436
	437	=head1 SYNOPSIS
	438
	439	use Math::BigInt;
	440	$i = Math::BigInt->new($string);
	441
	442	$i->bneg return BINT negation
	443	$i->babs return BINT absolute value
	444	$i->bcmp(BINT) return CODE compare numbers (undef,<0,=0,>0)
	445	$i->badd(BINT) return BINT addition
	446	$i->bsub(BINT) return BINT subtraction
	447	$i->bmul(BINT) return BINT multiplication
	448	$i->bdiv(BINT) return (BINT,BINT) division (quo,rem) just quo if scalar
	449	$i->bmod(BINT) return BINT modulus
	450	$i->bgcd(BINT) return BINT greatest common divisor
	451	$i->bnorm return BINT normalization
e16b8f49	452	$i->blsft(BINT) return BINT left shift
	453	$i->brsft(BINT) return (BINT,BINT) right shift (quo,rem) just quo if scalar
	454	$i->band(BINT) return BINT bit-wise and
	455	$i->bior(BINT) return BINT bit-wise inclusive or
	456	$i->bxor(BINT) return BINT bit-wise exclusive or
	457	$i->bnot return BINT bit-wise not
a5f75d66	458
	459	=head1 DESCRIPTION
	460
	461	All basic math operations are overloaded if you declare your big
	462	integers as
	463
	464	$i = new Math::BigInt '123 456 789 123 456 789';
	465
	466
	467	=over 2
	468
	469	=item Canonical notation
	470
	471	Big integer value are strings of the form C</^[+-]\d+$/> with leading
	472	zeros suppressed.
	473
	474	=item Input
	475
	476	Input values to these routines may be strings of the form
	477	C</^\s*[+-]?[\d\s]+$/>.
	478
	479	=item Output
	480
	481	Output values always always in canonical form
	482
	483	=back
	484
	485	Actual math is done in an internal format consisting of an array
	486	whose first element is the sign (/^[+-]$/) and whose remaining
	487	elements are base 100000 digits with the least significant digit first.
	488	The string 'NaN' is used to represent the result when input arguments
	489	are not numbers, as well as the result of dividing by zero.
	490
	491	=head1 EXAMPLES
	492
	493	'+0' canonical zero value
	494	' -123 123 123' canonical value '-123123123'
	495	'1 23 456 7890' canonical value '+1234567890'
	496
	497
b3ac6de7	498	=head1 Autocreating constants
	499
	500	After C<use Math::BigInt ':constant'> all the integer decimal constants
e3c7ef20	501	in the given scope are converted to C<Math::BigInt>. This conversion
b3ac6de7	502	happens at compile time.
	503
	504	In particular
	505
	506	perl -MMath::BigInt=:constant -e 'print 2**100'
	507
8dcee03e	508	print the integer value of C<2**100>. Note that without conversion of
8dcee03e	509	constants the expression 2**100 will be calculated as floating point number.
b3ac6de7	510
a5f75d66	511	=head1 BUGS
	512
	513	The current version of this module is a preliminary version of the
	514	real thing that is currently (as of perl5.002) under development.
	515
	516	=head1 AUTHOR
	517
	518	Mark Biggar, overloaded interface by Ilya Zakharevich.
	519
	520	=cut