4 # Arbitrary size rational math package
8 # Input values to these routines consist of strings of the form
9 # m|^\s*[+-]?[\d\s]+(/[\d\s]+)?$|.
11 # "+0/1" canonical zero value
12 # "3" canonical value "+3/1"
13 # " -123/123 123" canonical value "-1/1001"
14 # "123 456/7890" canonical value "+20576/1315"
15 # Output values always include a sign and no leading zeros or
17 # This package makes use of the bigint package.
18 # The string 'NaN' is used to represent the result when input arguments
19 # that are not numbers, as well as the result of dividing by zero and
20 # the sqrt of a negative number.
21 # Extreamly naive algorthims are used.
23 # Routines provided are:
25 # rneg(RAT) return RAT negation
26 # rabs(RAT) return RAT absolute value
27 # rcmp(RAT,RAT) return CODE compare numbers (undef,<0,=0,>0)
28 # radd(RAT,RAT) return RAT addition
29 # rsub(RAT,RAT) return RAT subtraction
30 # rmul(RAT,RAT) return RAT multiplication
31 # rdiv(RAT,RAT) return RAT division
32 # rmod(RAT) return (RAT,RAT) integer and fractional parts
33 # rnorm(RAT) return RAT normalization
34 # rsqrt(RAT, cycles) return RAT square root
36 # Convert a number to the canonical string form m|^[+-]\d+/\d+|.
37 sub main'rnorm { #(string) return rat_num
40 if (m#^([+-]?\d+)(/(\d*[1-9]0*))?$#) {
41 &norm($1, $3 ? $3 : '+1');
47 # Normalize by reducing to lowest terms
48 sub norm { #(bint, bint) return rat_num
49 local($num,$dom) = @_;
52 } elsif ($dom eq 'NaN') {
54 } elsif ($dom =~ /^[+-]?0+$/) {
57 local($gcd) = &'bgcd($num,$dom);
59 $num = &'bdiv($num,$gcd);
60 $dom = &'bdiv($dom,$gcd);
65 substr($dom,0,1) = '';
71 sub main'rneg { #(rat_num) return rat_num
72 local($_) = &'rnorm($_[0]);
73 tr/-+/+-/ if ($_ ne '+0/1');
78 sub main'rabs { #(rat_num) return $rat_num
79 local($_) = &'rnorm($_[0]);
80 substr($_,0,1) = '+' unless $_ eq 'NaN';
85 sub main'rmul { #(rat_num, rat_num) return rat_num
86 local($xn,$xd) = split('/',&'rnorm($_[0]));
87 local($yn,$yd) = split('/',&'rnorm($_[1]));
88 &norm(&'bmul($xn,$yn),&'bmul($xd,$yd));
92 sub main'rdiv { #(rat_num, rat_num) return rat_num
93 local($xn,$xd) = split('/',&'rnorm($_[0]));
94 local($yn,$yd) = split('/',&'rnorm($_[1]));
95 &norm(&'bmul($xn,$yd),&'bmul($xd,$yn));
99 sub main'radd { #(rat_num, rat_num) return rat_num
100 local($xn,$xd) = split('/',&'rnorm($_[0]));
101 local($yn,$yd) = split('/',&'rnorm($_[1]));
102 &norm(&'badd(&'bmul($xn,$yd),&'bmul($yn,$xd)),&'bmul($xd,$yd));
106 sub main'rsub { #(rat_num, rat_num) return rat_num
107 local($xn,$xd) = split('/',&'rnorm($_[0]));
108 local($yn,$yd) = split('/',&'rnorm($_[1]));
109 &norm(&'bsub(&'bmul($xn,$yd),&'bmul($yn,$xd)),&'bmul($xd,$yd));
113 sub main'rcmp { #(rat_num, rat_num) return cond_code
114 local($xn,$xd) = split('/',&'rnorm($_[0]));
115 local($yn,$yd) = split('/',&'rnorm($_[1]));
116 &bigint'cmp(&'bmul($xn,$yd),&'bmul($yn,$xd));
120 sub main'rmod { #(rat_num) return (rat_num,rat_num)
121 local($xn,$xd) = split('/',&'rnorm($_[0]));
122 local($i,$f) = &'bdiv($xn,$xd);
130 # square root by Newtons method.
131 # cycles specifies the number of iterations default: 5
132 sub main'rsqrt { #(fnum_str[, cycles]) return fnum_str
133 local($x, $scale) = (&'rnorm($_[0]), $_[1]);
136 } elsif ($x =~ /^-/) {
139 local($gscale, $guess) = (0, '+1/1');
140 $scale = 5 if (!$scale);
141 while ($gscale++ < $scale) {
142 $guess = &'rmul(&'radd($guess,&'rdiv($x,$guess)),"+1/2");
144 "$guess"; # quotes necessary due to perl bug