Welcome, guest | Sign In | My Account | Store | Cart

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 Scott David Daniels (author) 20 years, 9 months ago

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":

n,d = farey(probability,lim)
print "Odds are %d : %d" % (n,d-n)

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:

if v < 0:
...
elif v < 0.5:
n,d = farey((v-v+1)/v, lim) # lim is wrong: decide what you want
return (d,n)
elif v > 1:
intpart = floor(v)
n,d = farey(v-intpart)
return n+intpart*d, d
...

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. Connelly Barnes 16 years, 3 months ago

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. Dick Moores 15 years ago

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 rdmoores@gmail.com Peter Nau 13 years, 3 months ago

Nice work. Thank you: I plan to reuse this code to convert probability floats to fractional representation (similar to the "odds" comment and example).

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. jerfelix 9 years, 12 months ago

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).

from __future__ import division

def not_farey(v,lim):
d= min((abs(int(round(v*d))/d-v),d) for d in xrange(1,lim))
return int(round(v*d)),d Created by Scott David Daniels on Mon, 2 Apr 2001 (PSF)

Required Modules

• (none specified)