lib/bigrat.pl

   1 package bigrat;
   2 require "bigint.pl";
   3 #
   4 # This library is no longer being maintained, and is included for backward
   5 # compatibility with Perl 4 programs which may require it.
   6 #
   7 # In particular, this should not be used as an example of modern Perl
   8 # programming techniques.
   9 #
  10 # Arbitrary size rational math package
  11 #
  12 # by Mark Biggar
  13 #
  14 # Input values to these routines consist of strings of the form
  15 #   m|^\s*[+-]?[\d\s]+(/[\d\s]+)?$|.
  16 # Examples:
  17 #   "+0/1"                          canonical zero value
  18 #   "3"                             canonical value "+3/1"
  19 #   "   -123/123 123"               canonical value "-1/1001"
  20 #   "123 456/7890"                  canonical value "+20576/1315"
  21 # Output values always include a sign and no leading zeros or
  22 #   white space.
  23 # This package makes use of the bigint package.
  24 # The string 'NaN' is used to represent the result when input arguments
  25 #   that are not numbers, as well as the result of dividing by zero and
  26 #       the sqrt of a negative number.
  27 # Extreamly naive algorthims are used.
  28 #
  29 # Routines provided are:
  30 #
  31 #   rneg(RAT) return RAT                negation
  32 #   rabs(RAT) return RAT                absolute value
  33 #   rcmp(RAT,RAT) return CODE           compare numbers (undef,<0,=0,>0)
  34 #   radd(RAT,RAT) return RAT            addition
  35 #   rsub(RAT,RAT) return RAT            subtraction
  36 #   rmul(RAT,RAT) return RAT            multiplication
  37 #   rdiv(RAT,RAT) return RAT            division
  38 #   rmod(RAT) return (RAT,RAT)          integer and fractional parts
  39 #   rnorm(RAT) return RAT               normalization
  40 #   rsqrt(RAT, cycles) return RAT       square root
  41 \f
  42 # Convert a number to the canonical string form m|^[+-]\d+/\d+|.
  43 sub main'rnorm { #(string) return rat_num
  44     local($_) = @_;
  45     s/\s+//g;
  46     if (m#^([+-]?\d+)(/(\d*[1-9]0*))?$#) {
  47         &norm($1, $3 ? $3 : '+1');
  48     } else {
  49         'NaN';
  50     }
  51 }
  52
  53 # Normalize by reducing to lowest terms
  54 sub norm { #(bint, bint) return rat_num
  55     local($num,$dom) = @_;
  56     if ($num eq 'NaN') {
  57         'NaN';
  58     } elsif ($dom eq 'NaN') {
  59         'NaN';
  60     } elsif ($dom =~ /^[+-]?0+$/) {
  61         'NaN';
  62     } else {
  63         local($gcd) = &'bgcd($num,$dom);
  64         $gcd =~ s/^-/+/;
  65         if ($gcd ne '+1') {
  66             $num = &'bdiv($num,$gcd);
  67             $dom = &'bdiv($dom,$gcd);
  68         } else {
  69             $num = &'bnorm($num);
  70             $dom = &'bnorm($dom);
  71         }
  72         substr($dom,$[,1) = '';
  73         "$num/$dom";
  74     }
  75 }
  76
  77 # negation
  78 sub main'rneg { #(rat_num) return rat_num
  79     local($_) = &'rnorm(@_);
  80     tr/-+/+-/ if ($_ ne '+0/1');
  81     $_;
  82 }
  83
  84 # absolute value
  85 sub main'rabs { #(rat_num) return $rat_num
  86     local($_) = &'rnorm(@_);
  87     substr($_,$[,1) = '+' unless $_ eq 'NaN';
  88     $_;
  89 }
  90
  91 # multipication
  92 sub main'rmul { #(rat_num, rat_num) return rat_num
  93     local($xn,$xd) = split('/',&'rnorm($_[$[]));
  94     local($yn,$yd) = split('/',&'rnorm($_[$[+1]));
  95     &norm(&'bmul($xn,$yn),&'bmul($xd,$yd));
  96 }
  97
  98 # division
  99 sub main'rdiv { #(rat_num, rat_num) return rat_num
 100     local($xn,$xd) = split('/',&'rnorm($_[$[]));
 101     local($yn,$yd) = split('/',&'rnorm($_[$[+1]));
 102     &norm(&'bmul($xn,$yd),&'bmul($xd,$yn));
 103 }
 104 \f
 105 # addition
 106 sub main'radd { #(rat_num, rat_num) return rat_num
 107     local($xn,$xd) = split('/',&'rnorm($_[$[]));
 108     local($yn,$yd) = split('/',&'rnorm($_[$[+1]));
 109     &norm(&'badd(&'bmul($xn,$yd),&'bmul($yn,$xd)),&'bmul($xd,$yd));
 110 }
 111
 112 # subtraction
 113 sub main'rsub { #(rat_num, rat_num) return rat_num
 114     local($xn,$xd) = split('/',&'rnorm($_[$[]));
 115     local($yn,$yd) = split('/',&'rnorm($_[$[+1]));
 116     &norm(&'bsub(&'bmul($xn,$yd),&'bmul($yn,$xd)),&'bmul($xd,$yd));
 117 }
 118
 119 # comparison
 120 sub main'rcmp { #(rat_num, rat_num) return cond_code
 121     local($xn,$xd) = split('/',&'rnorm($_[$[]));
 122     local($yn,$yd) = split('/',&'rnorm($_[$[+1]));
 123     &bigint'cmp(&'bmul($xn,$yd),&'bmul($yn,$xd));
 124 }
 125
 126 # int and frac parts
 127 sub main'rmod { #(rat_num) return (rat_num,rat_num)
 128     local($xn,$xd) = split('/',&'rnorm(@_));
 129     local($i,$f) = &'bdiv($xn,$xd);
 130     if (wantarray) {
 131         ("$i/1", "$f/$xd");
 132     } else {
 133         "$i/1";
 134     }
 135 }
 136
 137 # square root by Newtons method.
 138 #   cycles specifies the number of iterations default: 5
 139 sub main'rsqrt { #(fnum_str[, cycles]) return fnum_str
 140     local($x, $scale) = (&'rnorm($_[$[]), $_[$[+1]);
 141     if ($x eq 'NaN') {
 142         'NaN';
 143     } elsif ($x =~ /^-/) {
 144         'NaN';
 145     } else {
 146         local($gscale, $guess) = (0, '+1/1');
 147         $scale = 5 if (!$scale);
 148         while ($gscale++ < $scale) {
 149             $guess = &'rmul(&'radd($guess,&'rdiv($x,$guess)),"+1/2");
 150         }
 151         "$guess";          # quotes necessary due to perl bug
 152     }
 153 }
 154
 155 1;