Rational Arithmetic Computing with rational numbers opens new vistas. For example: +∘÷ \ 10 ⍴ 1 1 2 1.5 1.66667 1.6 1.625 1.61538 1.61905 1.61765 1.61818 ⍕Q (+Q)∘(÷Q) \ 10 ⍴ qi 1 1 2 3r2 5r3 8r5 13r8 21r13 34r21 55r34 89r55 We presents a facility for extended-precision rational arithmetic: a representation for rational numbers, a function qi (“rational input”) to specify rational numbers conveniently, and a monadic operator Q for working with rational numbers, This text, the slides, and the code are available for download from the Dyalog ’20 webpage. An earlier version of this text appeared in the Dyalog APL Chat Forum post Rational Arithmetic [Hui 2020]. A more detailed description of the implementation can be found in Rational Arithmetic Implementation Notes [Hui 2020b] and in Rational Arithmetic — RNG [Hui 2020c]. The original design was changed, informed by the analysis in A Design Exercise [Hui 2020d]. The operator “rats”in the dfns workspace [Scholes 2014]
also implements rational arithmetic;
the current facility is more general and more capable.
There was a Conference Edition of Dyalog APL
distributed at the 2011 Dyalog User Conference [Hui 2011],
in which rational arithmetic was implemented as part of the interpreter.
The facility was never released as part of the official interpreter
due to difficulties with making exact and inexact numbers
play nice together.
A rational number is a scalar of an enclosed 3-element vector, whose items are:
The numerator and the denominator are relatively prime. The vectors of digits are required to have no leading 0s (except for rational 0 , represented as ⊂(,¨0 1),0). Therefore, the representation is unique: rational numbers p and q are equal if and only if p≡q . Because a rational number so represented is a scalar, “structural” (⍴ , ↑ , etc.) and “selection” (⌿ , [;] , ⌷ , etc.) functions just work on arrays of rational numbers. For example: Hilbert←{÷1+∘.+⍨⍳⍵} Hilbert 5 ⍝ an array of floating point numbers 1 0.5 0.333333 0.25 0.2 0.5 0.333333 0.25 0.2 0.166667 0.333333 0.25 0.2 0.166667 0.142857 0.25 0.2 0.166667 0.142857 0.125 0.2 0.166667 0.142857 0.125 0.111111 qi Hilbert 5 ⍝ an array of rational numbers ┌───────┬───────┬───────┬───────┬───────┐ │┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│ ││1│1│1│││1│2│1│││1│3│1│││1│4│1│││1│5│1││ │└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│ ├───────┼───────┼───────┼───────┼───────┤ │┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│ ││1│2│1│││1│3│1│││1│4│1│││1│5│1│││1│6│1││ │└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│ ├───────┼───────┼───────┼───────┼───────┤ │┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│ ││1│3│1│││1│4│1│││1│5│1│││1│6│1│││1│7│1││ │└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│ ├───────┼───────┼───────┼───────┼───────┤ │┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│ ││1│4│1│││1│5│1│││1│6│1│││1│7│1│││1│8│1││ │└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│ ├───────┼───────┼───────┼───────┼───────┤ │┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│┌─┬─┬─┐│ ││1│5│1│││1│6│1│││1│7│1│││1│8│1│││1│9│1││ │└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│└─┴─┴─┘│ └───────┴───────┴───────┴───────┴───────┘ ⍕Q qi Hilbert 5 ⍝ rational formatted output 1 1r2 1r3 1r4 1r5 1r2 1r3 1r4 1r5 1r6 1r3 1r4 1r5 1r6 1r7 1r4 1r5 1r6 1r7 1r8 1r5 1r6 1r7 1r8 1r9 The function qi (“rational input”) provides a convenient way to specify rational numbers. It applies independently to items of the argument; each item can be one of several types:
qi has an optional left argument of the comparison tolerance to use (if different from the current ⎕ct value) in rationalizing a floating point number. Examples: qi ¯2.75 ┌──────────┐ │┌───┬─┬──┐│ ││1 1│4│¯1││ │└───┴─┴──┘│ └──────────┘ ⍕Q qi ¯2.75 ¯11r4 qi ⊂1 2 3 4 5 6 7 8 ┌─────────────────────┐ │┌───────────────┬─┬─┐│ ││1 2 3 4 5 6 7 8│1│1││ │└───────────────┴─┴─┘│ └─────────────────────┘ ⍕Q qi ⊂1 2 3 4 5 6 7 8 12345678 qi ○1 ⍝ ⎕ct=1e¯14 ┌───────────────────────────────┐ │┌─────────────┬─────────────┬─┐│ ││5 4 1 9 3 5 1│1 7 2 5 0 3 3│1││ │└─────────────┴─────────────┴─┘│ └───────────────────────────────┘ 18 ⍕ 5419351÷1725033 3.141592653589815___ 0 qi ○1 ⍝ ⎕ct=0 ┌─────────────────────────────────────┐ │┌─────────────────┬───────────────┬─┐│ ││2 4 5 8 5 0 9 2 2│7 8 2 5 6 7 7 9│1││ │└─────────────────┴───────────────┴─┘│ └─────────────────────────────────────┘ 18 ⍕ 245850922÷78256779 3.141592653589793___ 18 ⍕ ○1 3.141592653589793___ t←qi ⊂'1.2345e6789' ⍴¨¨t ⍝ # decimal digits in numerator & denominator ┌─────────┐ │┌────┬─┬┐│ ││6790│1│││ │└────┴─┴┘│ └─────────┘ ¯10 ⍕Q t ⍝ exponential format 1.234500000e6789 f Q derives a function which works on rational numbers and produces a rational result except as noted otherwise. Scalar Functions
If f and g are dyadic functions in the above list, then f Q/ , f Q⌿ , f Q\ , f Q⍀ , (f Q).(g Q) , and ∘.(f Q) also work on rational arguments. The boolean functions < ≤ = ≠ ≥ > return a boolean result,
an ordinary 0 or 1 instead
of a rational 0 or 1 .
Likewise ×⍵ returns an
ordinary ¯1 , 0 , or 1 .
Non-Scalar Functions ⊤Q ⍵ applies to rational scalar positive integer ⍵ and produces its binary representation. (⊤ does not have a monadic definition. One possibility is ⊤⍵ ←→ 2⊤⍵ , that is, bits or binary representation. (Similarly, ⊥⍵ would be 2⊥⍵ , binary value.) This is in fact the definition used in A Formal Description of SYSTEM/360 [Falkoff et al., 1964, Table 1]. It was actually implemented in APL\360 for a time but was later removed due to space limitations [Hui 1992]. For example, ⊤Q qi 123 ←→ 1 1 1 1 0 1 1 . ⍺ ⊤Q ⍵ applies to rational scalar positive integer ⍵ and produces its base-⍺ representation, for an ordinary integer ⍺ . For example, 5 ⊤Q qi 1234567 ←→ 3 0 4 0 0 1 2 3 2 . ⍕Q ⍵ formats rational array ⍵ , with an r separating the numerator and the denominator (if not 1). ⍺ ⍕Q ⍵ formats rational array ⍵ in fixed-point (if ⍺≥0) or exponential (if ⍺<0) format. The left argument ⍺ must be a single integer: ⍺≥0 specifies the number of digits after the decimal point; ⍺<0 specifies the number of decimal digits in exponential notation. ⌹Q is analoguous to ⌹ for rational arguments.
As usual, a DOMAIN ERROR is signalled if ⍵
is non-singular (but not as usual, if the error is signalled ⍵ really is singular).
Computing with rational numbers opens new vistas.
+∘÷\ 16⍴1 1 2 1.5 1.66667 1.6 1.625 1.61538 1.61905 1.61765 1.61818 1.61798 1.61806 1.61803 1.61804 1.61803 1.61803 ⍕Q (+Q)∘(÷Q)\ 16 ⍴ qi 1 1 2 3r2 5r3 8r5 13r8 21r13 34r21 55r34 89r55 144r89 233r144 377r233 610r377 987r610 1597r987 How about that: our friends the Fibonacci numbers. A001203← 3 7 15 1 292 1 1 1 2 1 3 1 14 ⍝ OEIS A001203 30 ⍕Q ⍪ (+Q)∘(÷Q)\ qi A001203 3.000000000000000000000000000000 3.142857142857142857142857142857 3.141509433962264150943396226415 3.141592920353982300884955752212 3.141592653011902604072261494774 3.141592653921421044708715941593 3.141592653467436705520454785349 3.141592653618936623397500301411 3.141592653581077771204419306582 3.141592653591403978482542414219 3.141592653589389171543687321707 3.141592653589815383241943777307 3.141592653589792659375626945712 Square Root via Newton’s Method If W is a number whose square root is sought and ⍵ is an approximation of the root, Newton’s method improves the approximation by computing the average of ⍵ , and W divided by ⍵ . By using rational arithmetic, a sequence of rational approximations obtains. sqrt←{ ⍝ apply Newton’s method ⍺ times to find sqrt ⍵ q2←qi 2 W←⍵ f←{q2÷Q⍨⍵+Q W÷Q ⍵} ⍺ {1≥⍺:⍵ ⋄ (⍺-1)∇ ⍵,f ¯1↑⍵} qi 1 } ⍕Q t←8 sqrt qi 2 1 3r2 17r12 577r408 665857r470832 886731088897r627013566048 1572584048032918633353217r1111984844349868137938112 ... The following code displays the rational approximations to 40 decimal places, and the differences between 2 and the squares of the approximations. The differences illustrate that Newton’s method doubles the number of correct digits with each iteration. (40 ⍕Q ⍪t), ' ', ¯6 ⍕Q ⍪ (qi 2) -Q t ×Q t 1.0000000000000000000000000000000000000000 1.00000e0 1.5000000000000000000000000000000000000000 ¯2.50000e¯1 1.4166666666666666666666666666666666666667 ¯6.94444e¯3 1.4142156862745098039215686274509803921569 ¯6.00730e¯6 1.4142135623746899106262955788901349101166 ¯4.51095e¯12 1.4142135623730950488016896235025302436150 ¯2.54358e¯24 1.4142135623730950488016887242096980785697 ¯8.08728e¯49 1.4142135623730950488016887242096980785697 ¯8.17550e¯98 Hilbert’s matrix is known to be invertible, but is notorious for being nearly singular. A← qi Hilbert 7 ⍕Q A 1 1r2 1r3 1r4 1r5 1r6 1r7 1r2 1r3 1r4 1r5 1r6 1r7 1r8 1r3 1r4 1r5 1r6 1r7 1r8 1r9 1r4 1r5 1r6 1r7 1r8 1r9 1r10 1r5 1r6 1r7 1r8 1r9 1r10 1r11 1r6 1r7 1r8 1r9 1r10 1r11 1r12 1r7 1r8 1r9 1r10 1r11 1r12 1r13 B← ⌹Q A ⍕Q A (+Q).(×Q) B 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 Inverting the 7-by-7 Hilbert matrix as 64-bit floating point numbers: X←Hilbert 7 Y←⌹X ¯3 ⍕ X +.× Y 1.00E0 1.46E¯11 1.16E¯10 4.66E¯10 9.31E¯10 ¯9.31E¯10 0.00E0 ¯9.09E¯13 1.00E0 ¯3.49E¯10 2.33E¯9 9.31E¯10 ¯4.66E¯10 0.00E0 0.00E0 ¯7.28E¯12 1.00E0 9.31E¯10 9.31E¯10 ¯1.40E¯9 4.66E¯10 2.27E¯13 ¯2.18E¯11 ¯1.16E¯10 1.00E0 1.86E¯9 ¯2.33E¯9 2.33E¯10 4.55E¯13 ¯2.91E¯11 ¯5.82E¯11 4.66E¯10 1.00E0 ¯2.79E¯9 5.82E¯10 0.00E0 2.91E¯11 ¯3.49E¯10 9.31E¯10 ¯1.40E¯9 1.00E0 ¯4.66E¯10 ¯4.55E¯13 ¯1.46E¯11 ¯1.16E¯10 9.31E¯10 ¯1.86E¯9 0.00E0 1.00E0 How does it look for a matrix of random numbers? ⎕rl←6*5 ⍝ for reproducible random numbers X←¯10+?6 6⍴25 ¯3 ⍕ X +.× ⌹X 1.00E0 3.55E¯15 ¯2.84E¯14 ¯1.78E¯14 2.49E¯14 ¯2.22E¯15 5.55E¯17 1.00E0 2.84E¯14 2.13E¯14 ¯7.11E¯15 0.00E0 ¯6.94E¯17 3.55E¯15 1.00E0 ¯2.13E¯14 2.13E¯14 2.11E¯15 2.78E¯17 3.55E¯15 7.11E¯15 1.00E0 0.00E0 ¯8.88E¯16 ¯1.39E¯16 ¯1.07E¯14 ¯1.42E¯14 1.42E¯14 1.00E0 8.88E¯16 3.47E¯16 ¯9.77E¯15 ¯1.78E¯14 2.13E¯14 ¯2.31E¯14 1.00E0 A← qi X ⍕Q A (+Q).(×Q) ⌹Q A 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 Addition examples of using rational numbers can be found in
the APL Chat Forum posts modpower [Hui 2020e]
and Miller-Rabin Primality Testing [Hui 2020f].
A more detailed description of the rational arithmetic implementation can be found in Rational Arithmetic Implementation Notes [Hui 2020b] and in Rational Arithmetic — RNG [Hui 2020c]. The original design was changed, informed by the analysis in A Design Exercise [Hui 2020d]. There was a “Conference Edition” of Dyalog APL distributed at the 2011 Dyalog User Conference [Hui 2011], in which rational arithmetic was implemented as part of the interpreter. The facility was never released as part of the official interpreter due to difficulties with making exact and inexact numbers play nice together. The operator “rats”in the dfns workspace [Scholes 2014]
also implements rational arithmetic;
the current facility is more general and more capable.
We have shown that it is possible to compute rational approximations to non-rational functions. Another approach is to have floating point arithmetic to a specified number of digits (“big floats”). The following examples are from work in progress. MAXDIGITS specifies the number of significant floating point decimal digits. Like the rational arithmetic facility, there is a specification for the representation, a function fi (“floating point input”) to specify big floats conveniently, and a monadic operator F for working with big floats. The examples illustrate the computation of Euler’s number (*1) using a Taylor series. (The number of terms in the series, 60, is estimated by !⍣¯1 .) MAXDIGITS ← 80 fi 3.14159e¯20 ┌───────────────────┐ │┌───────────┬───┬─┐│ ││3 1 4 1 5 9│¯25│1││ │└───────────┴───┴─┘│ └───────────────────┘ ¯20 ⍕ +/ ÷ ! ⍳60 2.718281828459046____E0 ¯80 ⍕F +F/ ÷F !F fi ⍳60 2.71828182845904523536028747135266249775724709369995957 476e0 ' ',⎕d[A001113] ⍝ OEIS A001113 271828182845904523536028747135266249775724709369995957 4759 !⍣¯1 ⊢ 1e80 58.92 References
|