[math-fun] Turning a fraction into an integer and reverse: three efficient algorithms
What we want is three algorithms. Algorithm 1 inputs integer N and outputs rational p/q. Algorithm 2 ("ordinal") inputs p/q and outputs integer N (inverse of alg 1). Algorithm 3 inputs rational p/q, and effectively converts it to integer, increments integer, and converts it back to rational which it outputs ("successor"), except it hopefully does it faster than that. The Stern sequence doubling recurrences in A002487 show how alg 1 can be implemented in a number of arithmetic ops proportional to the number of bits in the integer. As indeed Dijkstra pointed out in both https://www.cs.utexas.edu/users/EWD/transcriptions/EWD05xx/EWD570.html and https://www.cs.utexas.edu/users/EWD/transcriptions/EWD05xx/EWD578.html To redo Dijkstra my way, given integer N we may compute the rational p/q where p=fusc(N) and q=fusc(N+1) by letting (p,q)=Alg1(N) where, recursively, Alg1(N){ if(N>0){ if(N is even){ (p,q)=Alg1(N/2); return( p, p+q ); }else{ (p,q)=Alg1((N-1)/2); return( p+q, q ); } } return( 0,1 ); } This can be rewritten as an iterative nonrecursive algorithm: Alg1(N){ p=0; q=1; for( the bits b of N from most-signif to least-signif do ){ if(b){ p = p+q; }else{ q = p+q; } } return( p,q ); } That in turn tells us how to write the inverse algorithm Alg2(p,q){ N=0; for(j=1; p>0; j=2*j){ if(q>p){ q = q-p; }else{ p = p-q; N += j; } } return(N); } which note also runs in a number of arithmetic ops proportional to the number of bits in N. A002487 says Moshe Newman pointed out how to implement Alg3 in a constant number of arithmetic ops: the successor of x=p/q is s=1/(2*floor(x) + 1 - x) = 1/(floor(x) + 1 - frac(x)). Does this really count as a constant number of ops -- how hard is the GCD you need to reduce the new fraction? Well, letting p=q*t+u where 0<=u<q so that p/q = t + u/q, evidently GCD(u,q)=1 since p/q was already in lowest terms, so then we see s=q/(t*q+q-u)=q/(p+q-2*u) is already in lowest terms unless u=0. If u=0 (which only happens when q=1) then s=1/(t+1) is lowest terms. Hence, no GCDing ever needs to be done. This yields the very simple rational-successor algorithm Alg3(p,q){ u = (p mod q); return( q, p+q-2*u ); } --Warren D. Smith http://RangeVoting.org
participants (1)
-
Warren D Smith