This converts a Numeric to a rational. The result is always in reduced form, but the proof, while possible, is subtle.
farey(math.pi,100) = (22,7)
Python, 22 lines
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22
def farey(v, lim): '''Named after James Farey, an English surveyor. No error checking on args -- lim = max denominator, results are (numerator, denominator), (1,0) is infinity ''' if v < 0: n,d = farey(-v, lim) return -n,d z = lim-lim # get 0 of right type for denominator lower, upper = (z,z+1), (z+1,z) while 1: mediant = (lower + upper), (lower+upper) if v * mediant > mediant: if lim < mediant: return upper lower = mediant elif v * mediant == mediant: if lim >= mediant: return mediant if lower < upper: return lower return upper else: if lim < mediant: return lower upper = mediant
notes on farey. Note the trickiness with "z" -- it is a zero of the same type as the argument lim. This allows you to use longs as the limit if need be. To print "odds":
This code is ideally suited for re-implementation in a lower level language (say C or assembly) if you have the need or desire for rational or "odds" output. Because this uses only multiplication and addition, it can play to hardware strengths.
If you are using this in an environment where you call it with a lot of values very near 0., 1., or 0.5 (or _very_ simple fractions), you may find it too slow. You may improve its performance in a "continued fraction" style by appending to the first if:
James Farey was an English surveyor who wrote a letter to the "Journal of Science" around the end of the eighteenth century. In that letter he observed that, while reading a privately published list of the decimal equivalents of fractions, he noticed the following. For any three consecutive fractions in least terms (say, A/B, C/D, E/F), the middle one (C/D) is equal to the (A+E)/(B+F).
I enjoy my own image of Mr. Farey sitting up late on a rainy English night, reading tables of decimal expansions of fractions by an oil lamp. Calculation has come a long way since his day, and I'm pleased to be able to benefit from his work.
Continued fractions. Alternatively, you can use continued fractions: http://mathworld.wolfram.com/ContinuedFraction.html . Note that in the continued fraction, you can use subtractions as well as additions. Thus math.pi == 3+1/(7+1/(16-.003405593314)) ~= 3+1/(7+1/16.) == 355./113.
farey() often gives incorrect result. For example, farey(0.36, 10) returns (1, 3), whereas (3, 8) is correct.
For another, farey(0.584115140346, 100) returns (7, 12) whereas (52, 89) is correct.
Dick Moores email@example.com
Nice work. Thank you: I plan to reuse this code to convert probability floats to fractional representation (similar to the "odds" comment and example).
I have comments on the comment about farey(...) accuracy:
The code appears to make an approximation that is increasingly accurate as one enlarges lim.
Note that farey(0.36, 100) returns (9, 25), which is the exact answer. (3,8) is "less" correct.
It's instructive to consider the following examples:
farey(0.584115140346, 100) returns (7, 12)
farey(0.584115140346, 1000) returns (125, 214)
farey(0.584115140346, 10000) returns (5479, 9380)
farey(0.584115140346, 100000) returns (37419, 64061)
In fact, farey(0.5842696, 100) returns (52, 89), which is "correct" -- to a desired degree of approximation.
I hope that's helpful.
I tend to agree with Dick Moores, above, as opposed to Peter Nau. The second parameter specifies the limit to the denominator, and you'd think that the integer ratio which best matches the first parameter and keeps the denominator below the second parameter would be the match.
Here's some code to address the accuracy issue. My code uses the 2nd parameter (lim) as a value that the denominator must be strictly less than (as opposed to less than OR equal to).