[KBK] says: I liked RS's fraction procedures in [Bag of algorithms], but I didn't like the comment that "precision must be specified." I saw it as a challenge to come up with something that would do: ---- % num2fract 0.125 1/8 % num2fract 0.44 11/25 % num2fract 0.3333333333 1/3 % num2fract 0.5714285714 4/7 ---- In short, I wanted to come as close as possible to converting ANY float to a fraction. This actually isn't as hard as might seem. The following code does it reasonably well (and is much improved [[2002-01-11]] over an earlier version that lived on the Wiki for over a year). ---- # Convert a floating point number to a mixed number, e.g., # 3.75 -> 3 3/4. The optional arg "maxint" is the largest # integer that will appear in the denominator. proc num2mixed { float { maxint 2147483647 } } { set numer [expr { int( $float ) }] set denom [expr { $float - $numer }] return [list $numer [num2fract $denom $maxint]] } # Convert a floating point number to a fraction, e.g., # 3.375 -> 27/8. The optional arg "maxint" is the largest # integer that will appear in the fraction, e.g., # [num2fract 3.14159 1000] -> 355/113 proc num2fract { float { maxint 2147483647 } } { return [join [quotient_rep $float $maxint] "/"] } [KBK] (11 January 2002) - I've replaced quotient_rep with an 'improved' version, which, alas, doesn't use the tricky [[catch]]. The old version is at [Tricky catch]. # Convert a floating point number to a pair of integers # that are the numerator and denominator of its fractional # representation. The optional arg "maxint" is the largest # integer that will appear in the fraction, e.g., # [quotient_rep 3.14159 1000] -> {355 113} proc quotient_rep { num { maxint 2147483647 } } { set i 0 set p_i 1 set q_i 0 set p_im1 0 set q_im1 1 while { 1 } { incr i set a_i [expr { int($num) }] set fract [expr { $num - $a_i }] if { 1.0 * $a_i * $p_i + $p_im1 > $maxint } { break } if { 1.0 * $a_i * $q_i + $q_im1 > $maxint } { break } foreach { p_i p_im1 } \ [list [expr { $a_i * $p_i + $p_im1 }] $p_i] break foreach { q_i q_im1 } \ [list [expr { $a_i * $q_i + $q_im1 }] $q_i] break if { abs( $fract * $maxint ) < 1 } { break } set num [expr { 1.0 / $fract }] } return [list $p_i $q_i] } '''Lars Hellström''': The two '''foreach''' commands above are unnecessary. An alternative method for doing the same thing is set p_i [expr {$p_im1 + $a_i*[set p_im1 $p_i]}] set q_i [expr {$q_im1 + $a_i*[set q_im1 $q_i]}] which works because '''expr''' looks up the value of e.g. '''p_im1''' for '''$p_im1''' ''before'' the '''set p_im1 $p_i''' is evaluated. The reader who wants to know the maths behind the above (regular continued fractions and how this algorithm computes them) can take a look at [http://abel.math.umu.se/~lars/misc/pstofont-excerpt.pdf], where the necessary formulae are derived from elementary mathematics. It is an excerpt from the source of a PS type 1 font linker I've written. (Type 1 fonts probably provide one of the few practical applications of algorithms such as the above, since the Type 1 charstring format [http://partners.adobe.com/asn/developer/pdfs/tn/T1_SPEC.PDF] offers no other way of encoding a non-integer than as the fraction of two integers.) [KBK] ''Exercise'' (for the mathematically inclined): Prove that the [[while]] loop terminates no later than ''$i == 41'' and find an input to [[quotient_rep]] for which the performance is that bad. [Fraction Math - Answer to exercise] [Arjen Markus] I have not yet analysed the above (not even the original one, though both are intriguing), but I do know a two-step approximation that does seem nice as well. I will post it here as soon as convenient. The only drawback: no guarantee that is "best", only that it is better than the trivial first step. [MS] Note that this algorithm does produce an approx p/q that is "best" approximator to num in the following sense: ''no other rational r/s with s<=q is closer''. But it is not "best" in the sense: ''no other rational involving integers <= $maxint is closer": consider the counterexample quotient_rep 3.1416305 500 --> 355/113 However, the fraction 377/120 is a better approximation with both numerator and denominator <=500. [KBK] If you're curious what's going on here, and where else you can go with it, check out the wonderful musings at [http://www.inwap.com/pdp10/hbaker/hakmem/cf.html]! [MS] Some of the claims in there seem (at this point, to me) to be wrong - look at [Notes on continued fractions] for counterexamples and notes. [[''This paragraph has been slightly edited.'']] ---- ''RS'': Good job, as far I've seen it under time pressure. But just to justify why I had the mandatory denominator: I wanted to emulate the behavior on fractions of e.g. inches or stock quotes, where denominator was always a power of two. I think 2-1/3" is not really idiomatic ;-) ''DKF'': Just fixed a minor documentation typo in ''num2fract''. Now that we've got code to convert arbitrary decimal expansions to fractions, we can put together some simple '''quotient arithmetic''' code. Note that it tries to give you the same kind of fractions out as you put in. ''Cannot handle negative numbers, decimal input or very large quotients.'' Take this and build on it! namespace eval quotientmath { # First, some code to handle "p q/r" and "p/q" input and output # formats. I'll use a list format internally. proc Qscan {fraction} { set bits [split $fraction " "] if {[llength $bits] < 1 || [llength $bits] > 2} { return -code error "\"$fraction\" is not a fraction" } elseif {[llength $bits] == 1} { set int 0 set kind 0 } elseif {![string is integer [lindex $bits 0]]} { return -code error "\"$fraction\" is not a fraction" } else { set int [lindex $bits 0] set fraction [lindex $bits 1] set kind 1 } set bits [split $fraction "/"] if { [llength $bits] != 2 || ![string is integer [lindex $bits 0]] || ![string is integer [lindex $bits 1]] } then { return -code error "\"$fraction\" is not a fraction" } set numer [lindex $bits 0] set denom [lindex $bits 1] list $kind [list [expr {$numer+$int*$denom}] $denom] } proc Qformat {kind fraction} { foreach {numer denom} $fraction {} if {$kind && $numer>$denom} { set int [expr {$numer/$denom}] set numer [expr {$numer%$denom}] return [format "%d %d/%d" $int $numer $denom] } else { return [format "%d/%d" $numer $denom] } } # Now the basic arithmetic operators. proc Qmult {a b} { foreach {an ad} $a {}; foreach {bn bd} $b {} list [expr {$an*$bn}] [expr {$ad*$bd}] } proc Qdiv {a b} { foreach {an ad} $a {}; foreach {bn bd} $b {} list [expr {$an*$bd}] [expr {$ad*$bn}] } proc Qadd {a b} { foreach {an ad} $a {}; foreach {bn bd} $b {} list [expr {$an*$bd+$bn*$ad}] [expr {$ad*$bd}] } proc Qsub {a b} { foreach {an ad} $a {}; foreach {bn bd} $b {} list [expr {$an*$bd-$bn*$ad}] [expr {$ad*$bd}] } # A normalisation function proc Qnorm {x} { foreach {p q} $x {} # Handle canonical forms if {$p == 0} {return [list 0 1]} if {$q == 0} {return [list 1 0]} # Calculate GCD while {1} { if {$q == 0} { # Termination case break } elseif {$p>$q} { # Swap p and q set t $p set p $q set q $t } set q [expr {$q%$p}] } # Divide both numerator and denominator by the GCD foreach {xn xd} $x {} list [expr {$xn/$p}] [expr {$xd/$p}] } # The public functions proc add {x y} { foreach {kx x} [Qscan $x] {} foreach {ky y} [Qscan $y] {} Qformat [expr {$kx||$ky}] [Qnorm [Qadd $x $y]] } proc sub {x y} { foreach {kx x} [Qscan $x] {} foreach {ky y} [Qscan $y] {} Qformat [expr {$kx||$ky}] [Qnorm [Qsub $x $y]] } proc mult {x y} { foreach {kx x} [Qscan $x] {} foreach {ky y} [Qscan $y] {} Qformat [expr {$kx||$ky}] [Qnorm [Qmult $x $y]] } proc div {x y} { foreach {kx x} [Qscan $x] {} foreach {ky y} [Qscan $y] {} Qformat [expr {$kx||$ky}] [Qnorm [Qdiv $x $y]] } namespace export add sub mult div } ---- [RS] has this tiny "rationalizer": proc dbl2frac {dbl {eps 0.001}} { for {set den 1} {$den<1024} {incr den} { set num [expr {round($dbl*$den)}] if {abs(double($num)/$den - $dbl) < $eps} break } list $num $den } % dbl2frac 0.333 1 3 % dbl2frac 1.23 16 13 % dbl2frac -0.42 -13 31 ---- [Category Mathematics]