README

   1 = Continued Fraction Library =
   2 Ben Lynn, September 2008
   3
   4 == Introduction ==
   5
   6 1/3 has no finite floating point representation.  How can a standard tool for
   7 working with real numbers be unable to handle a number a young child could
   8 understand?
   9
  10 Enter continued fractions. Finite continued fractions uniquely (if you forbid
  11 the last term from being 1) and exactly represent rationals.  Contrast this
  12 with 1/3 = 0.333... and 1 = 0.999...  They live in a realm free from the
  13 tyranny of base 2. Commonly encountered irrational numbers, both algebraic and
  14 transcendental, are readily represented with infinite continued fractions, and
  15 when we truncate them we have an approximation along with an error bound. In
  16 contrast, floating point numbers must be accompanied at all times with their
  17 error bounds.
  18
  19 Repeated digits in decimal representation occur for fractions, but why waste
  20 such an interesting phenomenon on something as mundane as a rational? With
  21 continued fractions, repeated terms occur only for square roots. Other infinite
  22 but patterned sequences represent other celebrated numbers, such as e and pi.
  23
  24 Let's say we use floating-point to compute some answer to 100 decimal places.
  25 What if we now want 200 decimal places? Even if we possess the first 100
  26 digits, and carefully preserved the program state, we must start from scratch
  27 and redo all our arithmetic with 200 digit precision.
  28
  29 If only we had used continued fractions, for then we could arbitrarily
  30 increase the precision of our answer without redoing any work. Furthermore,
  31 the error analysis becomes trivial.
  32
  33 == Disadvantages ==
  34
  35 We must pay a price for the advantages of continued fractions.  We
  36 have a representation from which we can instantly infer error bounds. A
  37 representation for which we may arbitrarily increase precision without redoing
  38 any work. A representation which in some sense is the closest we can hope to
  39 get to an exact number. Is it any wonder therefore that binary operations,
  40 whose output must conform to this wondrous representation, are clumsy and
  41 convoluted, and that deriving continued fraction expansions for many
  42 reals is time-consuming?
  43
  44 These drawbacks are almost certainly why continued fractions remain obscure.
  45 Yet surely there are some problems that scream for continued fractions?
  46 Imagine a real-time application where on a timer interrupt, all continued
  47 fractions simply output the last computed convergent.  Under heavy system load,
  48 approximations are coarser but the show goes on.
  49
  50 == Installation ==
  51
  52 I use the GMP library for multi-precision arithmetic. Type
  53
  54   $ make
  55
  56 to compile. Try "pi 2000" to generate 2000 digits of pi using a continued
  57 fraction. It takes longer than I had hoped, but it's faster than
  58
  59   $ echo "scale = 2000; 4*a(1)" | bc -l
  60
  61 There's no API documentation; see the source.
  62
  63 == Demand channels ==
  64
  65 Threads are the future, if not already the present. Once, single core CPUs grew
  66 faster at an astounding pace. But limits are being approached, and engineers
  67 are turning to multicore designs to surmount technical obstacles. Multicore
  68 systems are already commonplace, and as time passes, the number of cores per
  69 system will steadily march upward. Happily, this suits continued fractions.
  70
  71 Inspired by Doug McIlroy's "Squinting at Power Series" (search for it), we
  72 spawn at least one thread per continued fraction, and also per operation on
  73 continued fractions. The threads communicate via crude demand channels, each
  74 providing a few output terms upon each request.
  75
  76 This scheme allows us to quickly write code for generating continued
  77 fraction expansions or operating on them. For example, consider the code for
  78 generating 'e' = 2.71828...
  79
  80 ..............................................................................
  81 // e = [2; 1, 2, 1, 1, 4, 1, ...]
  82 static void *e_expansion(cf_t cf) {
  83   int even = 2;
  84   cf_put_int(cf, even);
  85
  86   while(cf_wait(cf)) {
  87     cf_put_int(cf, 1);
  88     cf_put_int(cf, even);
  89     even += 2;
  90     cf_put_int(cf, 1);
  91   }
  92
  93   return NULL;
  94 }
  95 ..............................................................................
  96
  97 The function +cf_put+ places a term on the output channel, and +cf_wait+ is our
  98 implementation of demand channels, a sort of cooperative multitasking.  Without
  99 it, not only would we have to destroy and clean up after the thread ourselves,
 100 but more seriously, the thread might consume vast amounts of resources on
 101 unwanted terms.  The +cf_wait+ functions instructs the thread to stay idle.
 102 Our threads call this function often, and if it returns zero, our threads
 103 clean themselves up and exit.
 104
 105 == Bugs and other issues ==
 106
 107 I intend to devote time to projects that operate on real numbers rather than
 108 study real numbers for its own sake. But perhaps then I'll find a perfect
 109 application for continued fractions, which will push me to fix deficiencies
 110 in this code.
 111
 112 - The API could use cosmetic surgery.
 113
 114 - I want to replace the condition variable with a semaphore to silence Helgrind
 115   (whose documentation states that condition variables almost unconditionally
 116   and invariably raise false alarms).
 117
 118 - The Makefile is fragile and annoying to maintain.
 119
 120 - The testing infrastructure needs upgrading.
 121
 122 - I'd like a split function as described by McIlroy, though for continued
 123   fractions a "tee"-like function is more apt: rather than spawning a thread
 124   for each term held back, we spawn and watch over two threads on invocation to
 125   service the two demand channels, storing all backlogged terms in one place.