From 116a1b2f5ed9c18ce1d48b7a1e353c0db756b3e3 Mon Sep 17 00:00:00 2001 From: Steve Peters Date: Wed, 18 Jul 2007 18:52:12 +0000 Subject: [PATCH] Upgrade to Math-BigRat-0.20 p4raw-id: //depot/perl@31629 --- MANIFEST | 2 + lib/Math/BigRat.pm | 160 +++++++++++++++++++++++++++++++++++++++-- lib/Math/BigRat/t/bigfltpm.inc | 17 +++++ lib/Math/BigRat/t/bigfltrt.t | 4 +- lib/Math/BigRat/t/biglog.t | 99 +++++++++++++++++++++++++ lib/Math/BigRat/t/bigroot.t | 69 ++++++++++++++++++ 6 files changed, 344 insertions(+), 7 deletions(-) create mode 100644 lib/Math/BigRat/t/biglog.t create mode 100644 lib/Math/BigRat/t/bigroot.t diff --git a/MANIFEST b/MANIFEST index e7b7c0a..74494ae 100644 --- a/MANIFEST +++ b/MANIFEST @@ -2072,10 +2072,12 @@ lib/Math/BigRat.pm Math::BigRat lib/Math/BigRat/t/big_ap.t Math::BigRat test lib/Math/BigRat/t/bigfltpm.inc Math::BigRat test lib/Math/BigRat/t/bigfltrt.t Math::BigRat test +lib/Math/BigRat/t/biglog.t Math::BigRat test lib/Math/BigRat/t/bigratpm.inc Math::BigRat test lib/Math/BigRat/t/bigratpm.t Math::BigRat test lib/Math/BigRat/t/bigrat.t Math::BigRat test lib/Math/BigRat/t/bigratup.t test under $Math::BigInt::upgrade +lib/Math/BigRat/t/bigroot.t Math::BigRat test lib/Math/BigRat/t/requirer.t see if require works properly lib/Math/BigRat/t/trap.t see if trap_nan and trap_inf work lib/Math/Complex.pm A Complex package diff --git a/lib/Math/BigRat.pm b/lib/Math/BigRat.pm index 7732c365..a2f0796 100644 --- a/lib/Math/BigRat.pm +++ b/lib/Math/BigRat.pm @@ -23,7 +23,7 @@ use vars qw($VERSION @ISA $upgrade $downgrade @ISA = qw(Math::BigFloat); -$VERSION = '0.19'; +$VERSION = '0.20'; use overload; # inherit overload from Math::BigFloat @@ -338,6 +338,12 @@ sub config # return (later set?) configuration data as hash ref my $class = shift || 'Math::BigRat'; + if (@_ == 1 && ref($_[0]) ne 'HASH') + { + my $cfg = $class->SUPER::config(); + return $cfg->{$_[0]}; + } + my $cfg = $class->SUPER::config(@_); # now we need only to override the ones that are different from our parent @@ -1011,6 +1017,130 @@ sub blog $x->_new_from_float( $x->_as_float()->blog(Math::BigFloat->new("$y"),@r) ); } +sub bexp + { + # set up parameters + my ($self,$x,$y,$a,$p,$r) = (ref($_[0]),@_); + + # objectify is costly, so avoid it + if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) + { + ($self,$x,$y,$a,$p,$r) = objectify(2,$class,@_); + } + + return $x->binf() if $x->{sign} eq '+inf'; + return $x->bzero() if $x->{sign} eq '-inf'; + + # we need to limit the accuracy to protect against overflow + my $fallback = 0; + my ($scale,@params); + ($x,@params) = $x->_find_round_parameters($a,$p,$r); + + # also takes care of the "error in _find_round_parameters?" case + return $x if $x->{sign} eq 'NaN'; + + # no rounding at all, so must use fallback + if (scalar @params == 0) + { + # simulate old behaviour + $params[0] = $self->div_scale(); # and round to it as accuracy + $params[1] = undef; # P = undef + $scale = $params[0]+4; # at least four more for proper round + $params[2] = $r; # round mode by caller or undef + $fallback = 1; # to clear a/p afterwards + } + else + { + # the 4 below is empirical, and there might be cases where it's not enough... + $scale = abs($params[0] || $params[1]) + 4; # take whatever is defined + } + + return $x->bone(@params) if $x->is_zero(); + + # See the comments in Math::BigFloat on how this algorithm works. + # Basically we calculate A and B (where B is faculty(N)) so that A/B = e + + my $x_org = $x->copy(); + if ($scale <= 75) + { + # set $x directly from a cached string form + $x->{_n} = $MBI->_new("90933395208605785401971970164779391644753259799242"); + $x->{_d} = $MBI->_new("33452526613163807108170062053440751665152000000000"); + $x->{sign} = '+'; + } + else + { + # compute A and B so that e = A / B. + + # After some terms we end up with this, so we use it as a starting point: + my $A = $MBI->_new("90933395208605785401971970164779391644753259799242"); + my $F = $MBI->_new(42); my $step = 42; + + # Compute how many steps we need to take to get $A and $B sufficiently big + my $steps = Math::BigFloat::_len_to_steps($scale - 4); +# print STDERR "# Doing $steps steps for ", $scale-4, " digits\n"; + while ($step++ <= $steps) + { + # calculate $a * $f + 1 + $A = $MBI->_mul($A, $F); + $A = $MBI->_inc($A); + # increment f + $F = $MBI->_inc($F); + } + # compute $B as factorial of $steps (this is faster than doing it manually) + my $B = $MBI->_fac($MBI->_new($steps)); + +# print "A ", $MBI->_str($A), "\nB ", $MBI->_str($B), "\n"; + + $x->{_n} = $A; + $x->{_d} = $B; + $x->{sign} = '+'; + } + + # $x contains now an estimate of e, with some surplus digits, so we can round + if (!$x_org->is_one()) + { + # raise $x to the wanted power and round it in one step: + $x->bpow($x_org, @params); + } + else + { + # else just round the already computed result + delete $x->{_a}; delete $x->{_p}; + # shortcut to not run through _find_round_parameters again + if (defined $params[0]) + { + $x->bround($params[0],$params[2]); # then round accordingly + } + else + { + $x->bfround($params[1],$params[2]); # then round accordingly + } + } + if ($fallback) + { + # clear a/p after round, since user did not request it + delete $x->{_a}; delete $x->{_p}; + } + + $x; + } + +sub bnok + { + # set up parameters + my ($self,$x,$y,@r) = (ref($_[0]),@_); + + # objectify is costly, so avoid it + if ((!ref($_[0])) || (ref($_[0]) ne ref($_[1]))) + { + ($self,$x,$y,@r) = objectify(2,$class,@_); + } + + # do it with floats + $x->_new_from_float( $x->_as_float()->bnok(Math::BigFloat->new("$y"),@r) ); + } + sub _float_from_part { my $x = shift; @@ -1602,10 +1732,6 @@ Calculates the factorial of $x. For instance: Works currently only for integers. -=head2 blog() - -Is not yet implemented. - =head2 bround()/round()/bfround() Are not yet implemented. @@ -1739,6 +1865,30 @@ Compute $x ** $y. Please see the documentation in L for further details. +=head2 bexp() + + $x->bexp($accuracy); # calculate e ** X + +Calculates two integers A and B so that A/B is equal to C, where C is +Euler's number. + +This method was added in v0.20 of Math::BigRat (May 2007). + +See also L. + +=head2 bnok() + + $x->bnok($y); # x over y (binomial coefficient n over k) + +Calculates the binomial coefficient n over k, also called the "choose" +function. The result is equivalent to: + + ( n ) n! + | - | = ------- + ( k ) k!(n-k)! + +This method was added in v0.20 of Math::BigRat (May 2007). + =head2 config() use Data::Dumper; diff --git a/lib/Math/BigRat/t/bigfltpm.inc b/lib/Math/BigRat/t/bigfltpm.inc index 45f48ac..9006afe 100644 --- a/lib/Math/BigRat/t/bigfltpm.inc +++ b/lib/Math/BigRat/t/bigfltpm.inc @@ -115,6 +115,8 @@ while () $try .= '$x->facmp($y);'; } elsif ($f eq "fpow") { $try .= '$x ** $y;'; + } elsif ($f eq "bnok") { + $try .= '$x->bnok($y);'; } elsif ($f eq "froot") { $try .= "$setup; \$x->froot(\$y);"; } elsif ($f eq "fadd") { @@ -397,6 +399,21 @@ abc:+0:NaN +27:+90:270 +1034:+804:415668 $div_scale = 40; +&bnok ++inf:10:inf +NaN:NaN:NaN +NaN:1:NaN +1:NaN:NaN +1:1:1 +# k > n +1:2:0 +2:3:0 +# k < 0 +1:-2:0 +# 7 over 3 = 35 +7:3:35 +7:6:1 +100:90:17310309456440 &flog 0::NaN -1::NaN diff --git a/lib/Math/BigRat/t/bigfltrt.t b/lib/Math/BigRat/t/bigfltrt.t index 3ad98db..ff911b9 100755 --- a/lib/Math/BigRat/t/bigfltrt.t +++ b/lib/Math/BigRat/t/bigfltrt.t @@ -29,7 +29,7 @@ BEGIN plan tests => 1; } -use Math::BigRat::Test; # test via this Subclass +use Math::BigRat::Test lib => 'Calc'; # test via this Subclass use vars qw ($class $try $x $y $f @args $ans $ans1 $ans1_str $setup $CL); $class = "Math::BigRat::Test"; @@ -37,5 +37,5 @@ $CL = "Math::BigInt::Calc"; ok (1,1); -# fails stil 185 tests +# fails still too many tests #require 'bigfltpm.inc'; # all tests here for sharing diff --git a/lib/Math/BigRat/t/biglog.t b/lib/Math/BigRat/t/biglog.t new file mode 100644 index 0000000..d201c41 --- /dev/null +++ b/lib/Math/BigRat/t/biglog.t @@ -0,0 +1,99 @@ +#!/usr/bin/perl -w + +# Test blog function (and bpow, since it uses blog), as well as bexp(). + +use Test::More; +use strict; + +BEGIN + { + $| = 1; + # to locate the testing files + my $location = $0; $location =~ s/biglog.t//i; + if ($ENV{PERL_CORE}) + { + # testing with the core distribution + @INC = qw(../lib); + } + unshift @INC, '../lib'; + if (-d 't') + { + chdir 't'; + require File::Spec; + unshift @INC, File::Spec->catdir(File::Spec->updir, $location); + } + else + { + unshift @INC, $location; + } + print "# INC = @INC\n"; + + plan tests => 14; + } + +use Math::BigRat; + +my $cl = "Math::BigRat"; + +############################################################################# +# test log($n) + +# does not work yet +#is ($cl->new(2)->blog(), '0', "blog(2)"); +#is ($cl->new(288)->blog(), '5',"blog(288)"); +#is ($cl->new(2000)->blog(), '7', "blog(2000)"); + +############################################################################# +# test exp($n) + +is ($cl->new(1)->bexp()->as_int(), '2', "bexp(1)"); +is ($cl->new(2)->bexp()->as_int(), '7',"bexp(2)"); +is ($cl->new(3)->bexp()->as_int(), '20', "bexp(3)"); + +# rounding not implemented yet +#is ($cl->new(3)->bexp(10), '20', "bexp(3,10)"); + +# $x < 0 => NaN +ok ($cl->new(-2)->blog(), 'NaN'); +ok ($cl->new(-1)->blog(), 'NaN'); +ok ($cl->new(-10)->blog(), 'NaN'); +ok ($cl->new(-2,2)->blog(), 'NaN'); + +############################################################################# +# test bexp() with cached results + +is ($cl->new(1)->bexp(), + '90933395208605785401971970164779391644753259799242' . '/' . + '33452526613163807108170062053440751665152000000000', + 'bexp(1)'); +#is ($cl->new(2)->bexp(40), $cl->new(1)->bexp(45)->bpow(2,40), 'bexp(2)'); + +#is ($cl->new("12.5")->bexp(61), $cl->new(1)->bexp(65)->bpow(12.5,61), 'bexp(12.5)'); + +############################################################################# +# test bexp() with big values (non-cached) + +#is ($cl->new(1)->bexp(100), +# '2.718281828459045235360287471352662497757247093699959574966967627724076630353547594571382178525166427', +# 'bexp(100)'); + +is ($cl->new("12.5")->bexp(91), $cl->new(1)->bexp(95)->bpow(12.5,91), + 'bexp(12.5) to 91 digits'); + +############################################################################# +# some integer results +is ($cl->new(2)->bpow(32)->blog(2), '32', "2 ** 32"); +is ($cl->new(3)->bpow(32)->blog(3), '32', "3 ** 32"); +is ($cl->new(2)->bpow(65)->blog(2), '65', "2 ** 65"); + +my $x = Math::BigInt->new( '777' ) ** 256; +my $base = Math::BigInt->new( '12345678901234' ); +is ($x->copy()->blog($base), 56, 'blog(777**256, 12345678901234)'); + +$x = Math::BigInt->new( '777' ) ** 777; +$base = Math::BigInt->new( '777' ); +is ($x->copy()->blog($base), 777, 'blog(777**777, 777)'); + +# all done +1; + diff --git a/lib/Math/BigRat/t/bigroot.t b/lib/Math/BigRat/t/bigroot.t new file mode 100644 index 0000000..41fee89 --- /dev/null +++ b/lib/Math/BigRat/t/bigroot.t @@ -0,0 +1,69 @@ +#!/usr/bin/perl -w + +# Test broot function (and bsqrt() function, since it is used by broot()). + +# It is too slow to be simple included in bigfltpm.inc, where it would get +# executed 3 times. + +# But it is better to test the numerical functionality, instead of not testing +# it at all. + +use Test::More; +use strict; + +BEGIN + { + $| = 1; + # to locate the testing files + my $location = $0; $location =~ s/bigroot.t//i; + if ($ENV{PERL_CORE}) + { + # testing with the core distribution + @INC = qw(../lib); + } + unshift @INC, '../lib'; + if (-d 't') + { + chdir 't'; + require File::Spec; + unshift @INC, File::Spec->catdir(File::Spec->updir, $location); + } + else + { + unshift @INC, $location; + } + print "# INC = @INC\n"; + + plan tests => 4 * 2; + } + +use Math::BigFloat; +use Math::BigInt; + +my $cl = "Math::BigFloat"; +my $c = "Math::BigInt"; + +# 2 ** 240 = +# 1766847064778384329583297500742918515827483896875618958121606201292619776 + +# takes way too long +#test_broot ('2','240', 8, undef, '1073741824'); +#test_broot ('2','240', 9, undef, '106528681.3099908308759836475139583940127'); +#test_broot ('2','120', 9, undef, '10321.27324073880096577298929482324664787'); +#test_broot ('2','120', 17, undef, '133.3268493632747279600707813049418888729'); + +test_broot ('2','120', 8, undef, '32768'); +test_broot ('2','60', 8, undef, '181.0193359837561662466161566988413540569'); +test_broot ('2','60', 9, undef, '101.5936673259647663841091609134277286651'); +test_broot ('2','60', 17, undef, '11.54672461623965153271017217302844672562'); + +sub test_broot + { + my ($x,$n,$y,$scale,$result) = @_; + + my $s = $scale || 'undef'; + is ($cl->new($x)->bpow($n)->broot($y,$scale),$result, "Try: $cl $x->bpow($n)->broot($y,$s) == $result"); + $result =~ s/\..*//; + is ($c->new($x)->bpow($n)->broot($y,$scale),$result, "Try: $c $x->bpow($n)->broot($y,$s) == $result"); + } + -- 2.7.4