4 # Arbitrary size rational math package
6 # Input values to these routines consist of strings of the form
7 # m|^\s*[+-]?[\d\s]+(/[\d\s]+)?$|.
9 # "+0/1" canonical zero value
10 # "3" canonical value "+3/1"
11 # " -123/123 123" canonical value "-1/1001"
12 # "123 456/7890" canonical value "+20576/1315"
13 # Output values always include a sign and no leading zeros or
15 # This package makes use of the bigint package.
16 # The string 'NaN' is used to represent the result when input arguments
17 # that are not numbers, as well as the result of dividing by zero and
18 # the sqrt of a negative number.
19 # Extreamly naive algorthims are used.
21 # Routines provided are:
23 # rneg(RAT) return RAT negation
24 # rabs(RAT) return RAT absolute value
25 # rcmp(RAT,RAT) return CODE compare numbers (undef,<0,=0,>0)
26 # radd(RAT,RAT) return RAT addition
27 # rsub(RAT,RAT) return RAT subtraction
28 # rmul(RAT,RAT) return RAT multiplication
29 # rdiv(RAT,RAT) return RAT division
30 # rmod(RAT) return (RAT,RAT) integer and fractional parts
31 # rnorm(RAT) return RAT normalization
32 # rsqrt(RAT, cycles) return RAT square root
34 # Convert a number to the canonical string form m|^[+-]\d+/\d+|.
35 sub main'rnorm { #(string) return rat_num
38 if (m#^([+-]?\d+)(/(\d*[1-9]0*))?$#) {
39 &norm($1, $3 ? $3 : '+1');
45 # Normalize by reducing to lowest terms
46 sub norm { #(bint, bint) return rat_num
47 local($num,$dom) = @_;
50 } elsif ($dom eq 'NaN') {
52 } elsif ($dom =~ /^[+-]?0+$/) {
55 local($gcd) = &'bgcd($num,$dom);
57 $num = &'bdiv($num,$gcd);
58 $dom = &'bdiv($dom,$gcd);
63 substr($dom,0,1) = '';
69 sub main'rneg { #(rat_num) return rat_num
70 local($_) = &'rnorm($_[0]);
71 tr/-+/+-/ if ($_ ne '+0/1');
76 sub main'rabs { #(rat_num) return $rat_num
77 local($_) = &'rnorm($_[0]);
83 sub main'rmul { #(rat_num, rat_num) return rat_num
84 local($xn,$xd) = split('/',&'rnorm($_[0]));
85 local($yn,$yd) = split('/',&'rnorm($_[1]));
86 &norm(&'bmul($xn,$yn),&'bmul($xd,$yd));
90 sub main'rdiv { #(rat_num, rat_num) return rat_num
91 local($xn,$xd) = split('/',&'rnorm($_[0]));
92 local($yn,$yd) = split('/',&'rnorm($_[1]));
93 &norm(&'bmul($xn,$yd),&'bmul($xd,$yn));
97 sub main'radd { #(rat_num, rat_num) return rat_num
98 local($xn,$xd) = split('/',&'rnorm($_[0]));
99 local($yn,$yd) = split('/',&'rnorm($_[1]));
100 &norm(&'badd(&'bmul($xn,$yd),&'bmul($yn,$xd)),&'bmul($xd,$yd));
104 sub main'rsub { #(rat_num, rat_num) return rat_num
105 local($xn,$xd) = split('/',&'rnorm($_[0]));
106 local($yn,$yd) = split('/',&'rnorm($_[1]));
107 &norm(&'bsub(&'bmul($xn,$yd),&'bmul($yn,$xd)),&'bmul($xd,$yd));
111 sub main'rcmp { #(rat_num, rat_num) return cond_code
112 local($xn,$xd) = split('/',&'rnorm($_[0]));
113 local($yn,$yd) = split('/',&'rnorm($_[1]));
114 &bigint'cmp(&'bmul($xn,$yd),&'bmul($yn,$xd));
118 sub main'rmod { #(rat_num) return (rat_num,rat_num)
119 local($xn,$xd) = split('/',&'rnorm($_[0]));
120 local($i,$f) = &'bdiv($xn,$xd);
128 # square root by Newtons method.
129 # cycles specifies the number of iterations default: 5
130 sub main'rsqrt { #(fnum_str[, cycles]) return fnum_str
131 local($x, $scale) = (&'rnorm($_[0]), $_[1]);
134 } elsif ($x =~ /^-/) {
137 local($gscale, $guess) = (0, '+1/1');
138 $scale = 5 if (!$scale);
139 while ($gscale++ < $scale) {
140 $guess = &'rmul(&'radd($guess,&'rdiv($x,$guess)),"+1/2");
142 "$guess"; # quotes necessary due to perl bug