Real numbers into continued fractions

A rational approximation to a real number is an approximation in the form of the ratio of two integers, numerator and denominator. The denominator should be positive, and for definiteness we require the numerator and denominator to be relatively prime, that is with greatest common divisor equal to 1.

The simplest and most powerful to find rational approximation for real numbers is to express them as continued fractions, with the convergents forming the approximations. An algorithm to do this is most easily described as a series of steps.

0. Let x be a real number for which the approximation is wanted.

1. Write b₀ = ⌊x⌋ and put x = b₀ + (x − ⌊x⌋) = b₀ + r₀.

2. 0 ≤ r₀ < 1, by definition. There are two cases:

   • If r₀ = 0 the process is complete. The rational approximation is exact.

   • If r₀ > 0, note that 1/r₀ > 1. Write 1/r₀ = ⌊1/r₀⌋ + (1/r₀ − ⌊1/r₀⌋) = b₁ + r₁, with b₁ ≥ 1 and 0 ≤ r₁ < 1. Then: $$ x = b_0 + \frac{1}{1/r_0} = b_0 + \cfrac{1}{b_1 + r_1}$$

3. Continuing in this way we produce a continued fraction expansion for the real number of the form:

$$ x = b_0 + \cfrac{1}{b_1 + \cfrac{1}{b_2 + \cfrac{1}{b_3 + \dotsb \vphantom{\cfrac{1}{b_4}}}}} $$

where the partial numerators are all equal to 1 and the partial denominators are all positive integers. The fraction terminates if at any stage a ‘excess’ term, r_n, becomes 0, preventing the process from continuing.

Continued fractions with all a_i = 1, i = 1, 2, … are usually called simple.

F <- c(1, 1, numeric(15))
for(i in 3:17) ## Fibonacci sequence by recurrence
  F[i] <- F[i-1] + F[i-2]  
F

 [1]    1    1    2    3    5    8   13   21   34   55   89  144  233  377
[15]  610  987 1597

vfractional((sqrt(5) + 1)/2, eps = 0, maxConv = 1:16)

 [1] 1        2        3/2      5/3      8/5      13/8     21/13   
 [8] 34/21    55/34    89/55    144/89   233/144  377/233  610/377 
[15] 987/610  1597/987

partial_denominators <- function(x, k = 10) {
  b <- rep(NA, k)
  r <- x
  for(i in 1:k) {
    b[i] <- floor(r)
    r <- r - b[i]
    if(isTRUE(all.equal(r, 0))) break
    r <- 1/r
  }
  structure(b, names = paste0("b", 1:k-1))
}

x <- c(pi = base::pi, e = exp(1), phi = (sqrt(5) + 1)/2, 
       structure(sqrt(1:9), names = paste0("sqrt(", 1:9, ")")))
tab <- x %>% sapply(partial_denominators) %>% t 
tab[is.na(tab)] <- ""
kable(tab, align = "r", caption = "Partial denominators")

Some elementary properties

The convergents for a rational approximation constructed in this way have some interesting elementary properties. Again, all results are taken from Chapter 1 of Khovanskii (1963).

• For the convergents so obtained, P_n/Q_n, the numerators and denominators will be relatively prime, that is, the fraction will be expressed in its lowest terms.

• The even numbered convergents, P_2n/Q_2n form an increasing series of lower bounds to the true value and the odd numbered ones, P_2n + 1/Q_2n + 1 form a decreasing series of upper bounds: $$ \frac{P_0}{Q_0} < \frac{P_2}{Q_2} < \cdots < \frac{P_{2k}}{Q_{2k}} < x < \frac{P_{2k+1}}{Q_{2k+1}} < \cdots <\frac{P_3}{Q_3} < \frac{P_1}{Q_1} $$ Hence the errors at any stage, x − P_n/Q_n are alternately positive and negative.

• The absolute error at any stage is bounded as follows: $$\left|x - \frac{P_n}{Q_n}\right| < \frac{1}{Q_{n-1}Q_n}, \qquad n = 1, 2, \dots$$

• From the recurrence relation, the denominators, Q_n, are non-decreasing, ultimately increasing monotonically without limit, unless the continued fraction terminates after a finite number of convergents. It follows that the continued fraction either terminates at the true value, or converges in the limit to it. In practice, convergence is rapid.

We can illustrate some of these properties by looking at the sequence of convergents the process generates for π. Notice how the errors alternate in sign and rapidly decrease in absolute value. For reference, the accurate value of π is 3.141592653589793116.

Given a fixed denominator, d, the optimal rational approximation to x is clearly obtained by rounding xd, that is using $\left\lfloor xd+\frac12\right\rfloor/d$. The denominators Q_n of the convergents usually correspond to points where there is a sudden increase in accuracy. To illustrate this we consider approximating $\sqrt 5$. The first few convergents are as follows:

(s5 <- vfractional(sqrt(5), eps = 0, maxConv = 1:7))

[1] 2          9/4        38/17      161/72     682/305    2889/1292 
[7] 12238/5473

d5 <- denominators(s5)
e5 <- abs(sqrt(5) - numerical(s5))

Now consider all rational approximations with denominators no larger:

d <- seq(max(d5))
n <- round(sqrt(5) * d)

To simplify the graph we may remove any which are not in their lowest terms.

gcd <- mapply(FUN = function(a, b) if(b == 0) a else Recall(b, a %% b),
              n, d)
nd <- cbind(n, d)/gcd
nd <- nd[!duplicated(nd), ]
e <- abs(sqrt(5) - nd[, 1]/nd[, 2])

To see the detail of what is happening, we need to use a log-log scale:

The blue points give the errors in the approximations for all denominators up to 5473; the red points and lines indicate the subset of them which are the continued fraction convergents.

Note that the convergent, P_n/Q_n, is the most accurate approximation for any with denominator less that or equal to Q_n, and is closer than at most one or two with denominators less than Q_n + 1.