SYN SYN
[p5sagit/p5-mst-13.2.git] / lib / Math / BigFloat.pm
index 92e7016..74a023e 100644 (file)
@@ -4,40 +4,39 @@ use Math::BigInt;
 
 use Exporter;  # just for use to be happy
 @ISA = (Exporter);
+$VERSION = '0.01';     # never had version before
 
-%OVERLOAD = ( 
-                               # Anonymous subroutines:
-'+'    =>      sub {new BigFloat &fadd},
-'-'    =>      sub {new BigFloat
+use overload
+'+'    =>      sub {new Math::BigFloat &fadd},
+'-'    =>      sub {new Math::BigFloat
                       $_[2]? fsub($_[1],${$_[0]}) : fsub(${$_[0]},$_[1])},
-'<=>'  =>      sub {new BigFloat
-                      $_[2]? fcmp($_[1],${$_[0]}) : fcmp(${$_[0]},$_[1])},
-'cmp'  =>      sub {new BigFloat
-                      $_[2]? ($_[1] cmp ${$_[0]}) : (${$_[0]} cmp $_[1])},
-'*'    =>      sub {new BigFloat &fmul},
-'/'    =>      sub {new BigFloat 
+'<=>'  =>      sub {$_[2]? fcmp($_[1],${$_[0]}) : fcmp(${$_[0]},$_[1])},
+'cmp'  =>      sub {$_[2]? ($_[1] cmp ${$_[0]}) : (${$_[0]} cmp $_[1])},
+'*'    =>      sub {new Math::BigFloat &fmul},
+'/'    =>      sub {new Math::BigFloat 
                       $_[2]? scalar fdiv($_[1],${$_[0]}) :
                         scalar fdiv(${$_[0]},$_[1])},
-'neg'  =>      sub {new BigFloat &fneg},
-'abs'  =>      sub {new BigFloat &fabs},
+'neg'  =>      sub {new Math::BigFloat &fneg},
+'abs'  =>      sub {new Math::BigFloat &fabs},
 
 qw(
 ""     stringify
 0+     numify)                 # Order of arguments unsignificant
-);
+;
 
 sub new {
-  my $foo = fnorm($_[1]);
-  panic("Not a number initialized to BigFloat") if $foo eq "NaN";
-  bless \$foo;
+  my ($class) = shift;
+  my ($foo) = fnorm(shift);
+  bless \$foo, $class;
 }
+
 sub numify { 0 + "${$_[0]}" }  # Not needed, additional overhead
                                # comparing to direct compilation based on
                                # stringify
 sub stringify {
     my $n = ${$_[0]};
 
-    $n =~ s/^\+//;
+    my $minus = ($n =~ s/^([+-])// && $1 eq '-');
     $n =~ s/E//;
 
     $n =~ s/([-+]\d+)$//;
@@ -52,27 +51,13 @@ sub stringify {
     } else {
        $n = '.' . ("0" x (abs($e) - $ln)) . $n;
     }
+    $n = "-$n" if $minus;
 
     # 1 while $n =~ s/(.*\d)(\d\d\d)/$1,$2/;
 
     return $n;
 }
 
-# Arbitrary length float math package
-#
-# by Mark Biggar
-#
-# number format
-#   canonical strings have the form /[+-]\d+E[+-]\d+/
-#   Input values can have inbedded whitespace
-# Error returns
-#   'NaN'           An input parameter was "Not a Number" or 
-#                       divide by zero or sqrt of negative number
-# Division is computed to 
-#   max($div_scale,length(dividend)+length(divisor)) 
-#   digits by default.
-# Also used for default sqrt scale
-
 $div_scale = 40;
 
 # Rounding modes one of 'even', 'odd', '+inf', '-inf', 'zero' or 'trunc'.
@@ -84,27 +69,13 @@ sub fneg; sub fabs; sub fcmp;
 sub fround; sub ffround;
 sub fnorm; sub fsqrt;
 
-#   bigfloat routines
-#
-#   fadd(NSTR, NSTR) return NSTR            addition
-#   fsub(NSTR, NSTR) return NSTR            subtraction
-#   fmul(NSTR, NSTR) return NSTR            multiplication
-#   fdiv(NSTR, NSTR[,SCALE]) returns NSTR   division to SCALE places
-#   fneg(NSTR) return NSTR                  negation
-#   fabs(NSTR) return NSTR                  absolute value
-#   fcmp(NSTR,NSTR) return CODE             compare undef,<0,=0,>0
-#   fround(NSTR, SCALE) return NSTR         round to SCALE digits
-#   ffround(NSTR, SCALE) return NSTR        round at SCALEth place
-#   fnorm(NSTR) return (NSTR)               normalize
-#   fsqrt(NSTR[, SCALE]) return NSTR        sqrt to SCALE places
-
-\f
 # Convert a number to canonical string form.
 #   Takes something that looks like a number and converts it to
 #   the form /^[+-]\d+E[+-]\d+$/.
 sub fnorm { #(string) return fnum_str
     local($_) = @_;
     s/\s+//g;                               # strip white space
+    no warnings;       # $4 and $5 below might legitimately be undefined
     if (/^([+-]?)(\d*)(\.(\d*))?([Ee]([+-]?\d+))?$/ && "$2$4" ne '') {
        &norm(($1 ? "$1$2$4" : "+$2$4"),(($4 ne '') ? $6-length($4) : $6));
     } else {
@@ -115,6 +86,7 @@ sub fnorm { #(string) return fnum_str
 # normalize number -- for internal use
 sub norm { #(mantissa, exponent) return fnum_str
     local($_, $exp) = @_;
+       $exp = 0 unless defined $exp;
     if ($_ eq 'NaN') {
        'NaN';
     } else {
@@ -154,7 +126,7 @@ sub fmul { #(fnum_str, fnum_str) return fnum_str
        &norm(Math::BigInt::bmul($xm,$ym),$xe+$ye);
     }
 }
-\f
+
 # addition
 sub fadd { #(fnum_str, fnum_str) return fnum_str
     local($x,$y) = (fnorm($_[$[]),fnorm($_[$[+1]));
@@ -188,11 +160,12 @@ sub fdiv #(fnum_str, fnum_str[,scale]) return fnum_str
        $scale = length($xm)-1 if (length($xm)-1 > $scale);
        $scale = length($ym)-1 if (length($ym)-1 > $scale);
        $scale = $scale + length($ym) - length($xm);
-       &norm(&round(Math::BigInt::bdiv($xm.('0' x $scale),$ym),$ym),
+       &norm(&round(Math::BigInt::bdiv($xm.('0' x $scale),$ym),
+                   Math::BigInt::babs($ym)),
            $xe-$ye-$scale);
     }
 }
-\f
+
 # round int $q based on fraction $r/$base using $rnd_mode
 sub round { #(int_str, int_str, int_str) return int_str
     local($q,$r,$base) = @_;
@@ -203,12 +176,14 @@ sub round { #(int_str, int_str, int_str) return int_str
     } else {
        local($cmp) = Math::BigInt::bcmp(Math::BigInt::bmul($r,'+2'),$base);
        if ( $cmp < 0 ||
-                ($cmp == 0 &&
-                 ( $rnd_mode eq 'zero'                             ||
+                ($cmp == 0 &&                                          (
+                  ($rnd_mode eq 'zero'                            ) ||
                   ($rnd_mode eq '-inf' && (substr($q,$[,1) eq '+')) ||
                   ($rnd_mode eq '+inf' && (substr($q,$[,1) eq '-')) ||
-                  ($rnd_mode eq 'even' && $q =~ /[24680]$/)        ||
-                  ($rnd_mode eq 'odd'  && $q =~ /[13579]$/)        )) ) {
+                  ($rnd_mode eq 'even' && $q =~ /[13579]$/        ) ||
+                  ($rnd_mode eq 'odd'  && $q =~ /[24680]$/        )    )
+                 ) 
+               ) {
            $q;                     # round down
        } else {
            Math::BigInt::badd($q, ((substr($q,$[,1) eq '-') ? '-1' : '+1'));
@@ -233,7 +208,7 @@ sub fround { #(fnum_str, scale) return fnum_str
        }
     }
 }
-\f
+
 # round $x at the 10 to the $scale digit place
 sub ffround { #(fnum_str, scale) return fnum_str
     local($x,$scale) = (fnorm($_[$[]),$_[$[+1]);
@@ -248,7 +223,11 @@ sub ffround { #(fnum_str, scale) return fnum_str
            if ($xe < 1) {
                '+0E+0';
            } elsif ($xe == 1) {
-               &norm(&round('+0',"+0".substr($xm,$[+1,1),"+10"), $scale);
+               # The first substr preserves the sign, passing a non-
+               # normalized "-0" to &round when rounding -0.006 (for
+               # example), purely so &round won't lose the sign.
+               &norm(&round(substr($xm,$[,1).'0',
+                     "+0".substr($xm,$[+1,1),"+10"), $scale);
            } else {
                &norm(&round(substr($xm,$[,$xe),
                      "+0".substr($xm,$[+$xe,1),"+10"), $scale);
@@ -265,15 +244,24 @@ sub fcmp #(fnum_str, fnum_str) return cond_code
     if ($x eq "NaN" || $y eq "NaN") {
        undef;
     } else {
-       ord($y) <=> ord($x)
-       ||
-       (  local($xm,$xe,$ym,$ye) = split('E', $x."E$y"),
-            (($xe <=> $ye) * (substr($x,$[,1).'1')
-             || Math::BigInt::cmp($xm,$ym))
-       );
+       local($xm,$xe,$ym,$ye) = split('E', $x."E$y");
+       if ($xm eq '+0' || $ym eq '+0') {
+           return $xm <=> $ym;
+       }
+       if ( $xe < $ye )        # adjust the exponents to be equal
+       {
+               $ym .= '0' x ($ye - $xe);
+               $ye = $xe;
+       }
+       elsif ( $ye < $xe )     # same here
+       {
+               $xm .= '0' x ($xe - $ye);
+               $xe = $ye;
+       }
+       return Math::BigInt::cmp($xm,$ym);
     }
 }
-\f
+
 # square root by Newtons method.
 sub fsqrt { #(fnum_str[, scale]) return fnum_str
     local($x, $scale) = (fnorm($_[$[]), $_[$[+1]);
@@ -290,8 +278,91 @@ sub fsqrt { #(fnum_str[, scale]) return fnum_str
            $guess = fmul(fadd($guess,fdiv($x,$guess,$gs*2)),".5");
            $gs *= 2;
        }
-       new BigFloat &fround($guess, $scale);
+       new Math::BigFloat &fround($guess, $scale);
     }
 }
 
 1;
+__END__
+
+=head1 NAME
+
+Math::BigFloat - Arbitrary length float math package
+
+=head1 SYNOPSIS
+
+  use Math::BigFloat;
+  $f = Math::BigFloat->new($string);
+
+  $f->fadd(NSTR) return NSTR            addition
+  $f->fsub(NSTR) return NSTR            subtraction
+  $f->fmul(NSTR) return NSTR            multiplication
+  $f->fdiv(NSTR[,SCALE]) returns NSTR   division to SCALE places
+  $f->fneg() return NSTR                negation
+  $f->fabs() return NSTR                absolute value
+  $f->fcmp(NSTR) return CODE            compare undef,<0,=0,>0
+  $f->fround(SCALE) return NSTR         round to SCALE digits
+  $f->ffround(SCALE) return NSTR        round at SCALEth place
+  $f->fnorm() return (NSTR)             normalize
+  $f->fsqrt([SCALE]) return NSTR        sqrt to SCALE places
+
+=head1 DESCRIPTION
+
+All basic math operations are overloaded if you declare your big
+floats as
+
+    $float = new Math::BigFloat "2.123123123123123123123123123123123";
+
+=over 2
+
+=item number format
+
+canonical strings have the form /[+-]\d+E[+-]\d+/ .  Input values can
+have embedded whitespace.
+
+=item Error returns 'NaN'
+
+An input parameter was "Not a Number" or divide by zero or sqrt of
+negative number.
+
+=item Division is computed to 
+
+C<max($Math::BigFloat::div_scale,length(dividend)+length(divisor))>
+digits by default.
+Also used for default sqrt scale.
+
+=item Rounding is performed
+
+according to the value of
+C<$Math::BigFloat::rnd_mode>:
+
+  trunc     truncate the value
+  zero      round towards 0
+  +inf      round towards +infinity (round up)
+  -inf      round towards -infinity (round down)
+  even      round to the nearest, .5 to the even digit
+  odd       round to the nearest, .5 to the odd digit
+
+The default is C<even> rounding.
+
+=back
+
+=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.
+
+The printf subroutine does not use the value of
+C<$Math::BigFloat::rnd_mode> when rounding values for printing.
+Consequently, the way to print rounded values is
+to specify the number of digits both as an
+argument to C<ffround> and in the C<%f> printf string,
+as follows:
+
+  printf "%.3f\n", $bigfloat->ffround(-3);
+
+=head1 AUTHOR
+
+Mark Biggar
+
+=cut