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 However, we must pay a price for the advantages of continued fractions.  We
  34 have a representation from which we can instantly infer error bounds. A
  35 representation for which we may arbitrarily increase precision without redoing
  36 any work. A representation which in some sense is the closest we can hope to
  37 get to an exact number. Is it any wonder therefore that binary operations,
  38 whose output must conform to this wondrous representation, are clumsy and
  39 convoluted, and that deriving continued fraction expansions for many
  40 reals is time-consuming?
  41
  42 These drawbacks are almost certainly why continued fractions remain obscure.
  43 Yet surely there are some problems that scream for continued fractions?
  44 Imagine a real-time application where on a timer interrupt, all continued
  45 fractions simply output the last computed convergent.  Under heavy system load,
  46 approximations are coarser but the show goes on.
  47
  48 == Installation ==
  49
  50 I use the GMP library for multi-precision arithmetic. Type
  51
  52  $ make
  53
  54 to compile. The "pi" binary generates 5000 digits of pi using a continued
  55 fraction. It takes longer than I had hoped.
  56
  57 There's no API documentation yet; see the source.
  58
  59 == Demand channels ==
  60
  61 Threads are the future, if not already the present. Once, single core CPUs grew
  62 faster at an astounding pace. But limits are being approached, and engineers
  63 are turning to multicore designs to surmount technical obstacles. Multicore
  64 systems are already commonplace, and as time passes, the number of cores per
  65 system will steadily march upward. Happily, this suits continued fractions.
  66
  67 Inspired by Doug McIlroy's "Squinting at Power Series" (search for it), we
  68 spawn at least one thread per continued fraction, and also per operation on
  69 continued fractions. The threads communicate via crude demand channels, each
  70 providing a few output terms upon each request.
  71
  72 This scheme allows us to quickly write code for generating continued
  73 fraction expansions or operating on them. For example, consider the code for
  74 generating 'e' = 2.71828...
  75
  76 ..............................................................................
  77 // e = [2; 1, 2, 1, 1, 4, 1, ...]
  78 static void *e_expansion(cf_t cf) {
  79   int even = 2;
  80   cf_put_int(cf, even);
  81
  82   while(cf_wait(cf)) {
  83     cf_put_int(cf, 1);
  84     cf_put_int(cf, even);
  85     even += 2;
  86     cf_put_int(cf, 1);
  87   }
  88
  89   return NULL;
  90 }
  91 ..............................................................................
  92
  93 The function +cf_put+ places a term on the output channel, and +cf_wait+ is our
  94 implementation of demand channels, a sort of cooperative multitasking.  Without
  95 it, not only would we have to destroy and clean up after the thread ourselves,
  96 but more seriously, the thread might consume vast amounts of resources on
  97 unwanted terms.  The +cf_wait+ functions instructs the thread to stay idle.
  98 Our threads call this function often, and if it returns zero, our threads
  99 clean themselves up and exit.
 100
 101 == Bugs and other issues ==
 102
 103 I intend to devote time to projects that operate on real numbers rather than
 104 study real numbers for its own sake. But perhaps then I'll find a perfect
 105 application for continued fractions, which will push me to fix deficiencies
 106 in this code.
 107
 108 - The API could use cosmetic surgery.
 109
 110 - I want to replace the condition variable with a semaphore to silence Helgrind
 111   (whose documentation states that condition variables almost unconditionally
 112   and invariably raise false alarms).
 113
 114 - The Makefile is fragile and annoying to maintain.
 115
 116 - The testing infrastructure needs upgrading.
 117
 118 - I'd like a split function as described by McIlroy, though for continued
 119   fractions a "tee"-like function is more apt: rather than spawning a thread
 120   for each term held back, we spawn and watch over two threads on invocation to
 121   service the two demand channels, storing all backlogged terms in one thread.