I've been hacking a new version of the Common Lisp & Maxima "rationalize" function which produces the worst rational approximation to a complex number that falls within epsilon of that complex number. In the case of Maxima, the epsilon is set by the global variable "ratepsilon". The traditional Euclidean algorithm for this process focuses on the number itself (which I call c), and only later worries about the epsilon and the radius of convergence. My algorithm starts with the _circle_ whose center is the given complex number c and whose real radius is r; we can denote this circle by (c,r). Since translation and inversion in the complex plane both take circles to circles (more accurately, "generalized circles" to "generalized circles"), we can perform a Euclidean procedure on the _circle_ (c,r) itself. At each stage we keep track of the transformations (translations & inversions), so that we can recover the original circle, if necessary. We first look at the given circle (c,r). If there is a (Gaussian) integer point within this circle, then we are done, because such a point is rational, so this is our answer. If we compute cr=round(c) (= round(realpart(c))+i*round(imagpart(c))), then this number cr is (one of) the closest integer points to c. If cr is "inside" the circle (c,r), then we are done. Now if r>sqrt(2)/2, then cr _has_ to be inside the circle (c,r). So we now consider the case where cr is not inside (c,r). After translating the origin to cr, we know that this new center is within sqrt(2)/2 from the origin -- i.e., inside the unit circle. We then invert the given circle in the unit circle, which will put the center of that new circle outside the unit circle. Furthermore, the radius of that new circle will grow bigger. We can now compute the closest integral point to the center of that new circle, and so forth. At each of these Euclidean steps, the circle radius grows larger, and must eventually enclose an integral point, at which time we invert the sequence of transformations to compute the final rational answer. This answer will be rational, because the integral point was transformed by inversions and translations by (Gaussian) integers. This discussion is not a _proof_, but shows how the algorithm was conceived. Here's the algorithm in Maxima: euclid(zc,m):= block([q,iqmat], /* zc is input circle in Schwerdtfeger's representation */ /* m accumulates transforms; m starts out as the 2x2 identity matrix */ /* Original circle is always transpose(m^^-1).zc.conjugate(m^^-1) */ q:cround(center(zc)), if is(inside(q,zc)) then m.matrix([q],[1]) else (iqmat:matrix([1,q],[0,1]), euclid(transpose(iqmat.imat).zc.conjugate(iqmat.imat),m.iqmat.imat))); This produces is a column vector matrix([num],[denom]); the result is num/denom. euclid(gcircle(%pi,1.0e-7),ident(2)) produces [ - 122453 ] [ ] [ - 38978 ] = 122453/38978, which is within 10^-7 of pi. (%i13) abs((-122453)/(-38978)-%pi),numer; (%o13) 3.9724384226502707E-8 We use Hans Schwerdtfeger's representation for generalized circles as 2x2 Hermitian matrices. See this Wikipedia article for more details. http://en.wikipedia.org/wiki/Generalised_circle /* We test for "inside" by substituting the point into the equation. */ inside(w,c):=is(equation(c,w)<=0); Although very elegant, there is a problem with Schwerdtfeger's representation for circles: the squared radius of the circle is computed by the difference of two large quantities, so when we approach the limit of machine precision, we may inadvertently end up with a negative number, and hence an _imaginary_ radius! Perhaps this problem can be ameliorated with a slightly different representation.