5. Programs corresponding to recurrence relations

Theorem 8 mentions successive states connected by a recurrence relation. The meaning of this theorem is twofold: it can be used to prove assertions about a given program, but also —and this, I think, is more important— it suggests to us, when faced with the task of making a program, the use of a while-clause in the case of a problem that in its mathematical formulation presents itself as the evaluation of a recurrence relation. We are going to illustrate this by a number of examples.

Consider the sequence of pairs ai, ci given by

for i = 0 a0 = 1
c0 = 1 – b, with 0 < b < 2 (i.e. abs(c0) < 1)

for i > 0 ai = (1 + ci-1) * ai-1
ci = (ci-1)2     .

Then lim ai = 1/b     .
i →∞

Exercise. Prove the last formula. (This has nothing to do with programming, it is secondary school algebra. The clue of a proof can be found in the relation

1 1 + ci-1
-------   =   -------             .)
1 – ci-1 1 - ci

It is requested to use this recurrence relation to approximate the value of 1/b; obviously we cannot compute infinitely many elements of the sequence a0, a1, a2, … but we can accept ak as a sufficiently close (how close?) approximation of 1/b when ck is less (in absolute value) than a given, small, positive tolerance named "eps". (This example is of historical interest; it has been taken from the subroutine library for EDSAC 1, the world's first stored program controlled automatic computer. The order code of this computer did not comprise a divide instruction and one of the methods used with this computer to compute quotients was based on the above recurrence relation.)

Theorem 8 talks about "a part, s, of the state space" and the loop

while B(s) do s := S(s)

asserts that after the initial state s0 , the states si after the i-th execution of the repeatable statement will satisfy

si = S(si-1)

Our recurrence relations (2) are exactly of the form (3) if we identify the state si with the value pair ai, ci. That is, to span the part s of the state space we have to introduce two variables, for the purpose of this discussion called A and C, and we shall denote their values after the i-th execution of the repeatable statement Ai and Ci respectively. We associate the state si (as given by the values Ai and Ci ) with the value pair ai, ci by the relations

Ai = ai
Ci = ci

(Remember: on the left-hand sides the subscript "i" means "the value of the variable after the i-th execution of the repeatable statement", on the right-hand side the subscript "i" refers to the recurrent sequences as given by (1) and (2). It would have been usual to call the two variables "a" and "c" instead of "A" and "C", i.e. not to distinguish between the quantities defined in the mathematical formulation on the one hand and the associated variables in the program on the other hand. As this association is the very subject of this discussion, it would have been fatal not to distinguish between them.)

Within the scope of a declaration "real A, C" it is now a straightforward task to write the piece of program:

A := 1; C := 1 – b;
while abs(C ) ≥ eps do
begin A := (1 + C ) * A;
         C := C * C
end     .

The first line has to create the proper state s0 and does so in accordance with (4) and (1), the repeatable statement has the form, symbolically denoted by "s := S(s)" —see the Note below— in accordance with (4) and (2), and the condition guarantees that after termination

(Ak = ) A = ak

will hold with the proper value of k.

Exercise. Prove that the loop terminates.

Note. The symbolic assignment "s := S(s)" has the form of two assignments

A := (1 + C ) * A;
C := C * C     .

With the initial condition A = ai-1 , C = ci-1 the first assignment is equivalent to

A := (1 + ci-1) * ai-1

and after the first, but before the second assignment we have —on account of (2)—

A = ai, C = ci-1    .

We have the complete pair A = ai, C = ci only after the second assignment. Thanks to the explicit occurrence of the subscripts, the order of the two relations comprising (2) is immaterial, this in contrast to the two assignment statements composing the repeatable statement, whose order is vital.

Exercise. In the same EDSAC 1 subroutine library the next scheme is used. Consider the sequence of pairs ai, ci, given by

for i = 0 a0 = b
c0 = 1 – b, with 0 < b < 2 (i.e. abs(c0) < 1)

for i > 0 ai = (1 + .5 * ci-1) * ai-1
ci = (ci-1)2 * (.75 + .25 * ci-1)     .

Then lim ai = b.5     .

Prove the last formula and make a program using it for the approximation of the square root. The clue of a proof can be found in the relation

(1 – ci-1)-.5 = (1 + .5 * ci-1) * (1 - ci)-.5     .

Prove also the termination of the repetition in the program made.

Exercise. In the same EDSAC 1 subroutine library the next scheme is used. Consider the sequence of triples inci, si, xi, given by

for i = 0 inc0 = log 2
s0 = 0
x0 = arg               (with 1 ≤ arg < 2)
for i > 0
        for (xi-1)2 < 2 inci = .5 * inci-1
si = si-1
xi = (xi-1)2     .

       for (xi-1)2 ≥ 2 inci = .5 * inci-1
si = si-1 + .5 * inci-1
xi = .5 * (xi-1)2

Then lim si = log(arg)     .

Prove this relation and make a program using it to approximate the logarithm of a value arg in the interval stated. (In this program "log 2" may be regarded as a known constant which, in fact, determines the base of the logarithm.) The clue of the proof can be found in the relation

log(arg) = si + inci * log(xi ) / log 2     .

Our next example is very simple; it is so traditional that we could call it standard. (No self-respecting programming course omits it, it is often the very first example of a loop; Peter Naur uses it in his article "Proof of algorithms by general snapshots", 1966, BIT, 6, pp 310-316.)

Given a sequence of values

a[1], a[2], a[3], … , a[N]                  (with N ≥ 1)

and a variable called "max." Make a piece of program assigning to the variable named "max" the maximum value occurring in the sequence. (As N ≥ 1, the sequence is not empty and therefore the task makes sense; it is not required that any two values in the sequence differ from each other, the maximum value sought may occur more than once in the sequence.) If he welcomes the experience the reader is invited to try to make this piece of program himself before reading on.

How do we define the maximum value occurring in a sequence of length N for general N ≥ 1? If we call "maximumk" the maximum value occurring among the first k elements a[1], … , a[k], then

1) the answer sought is maximumN

2) the values maximumk are given

for k = 1 by the base: maximum1 = a[1]
appealing to the knowledge that the maximum element in a sequence of length 1 must be the only element in the sequence

for k > 1 by the recurrence relation:

maximumk = MAX(maximumk-1, a[k])

assuming the knowledge of the function MAX of two arguments.


The recurrence relation (6) presents us with an additional difficulty because it is not of the form

si = S(si-1)

because —via "a[k]"— the value k occurs on the right-hand side not exclusively in the subscript "k-1". To overcome this we use a trick that might be called a method. If we call nk the k-th natural number, then n= k; the numbers nk satisfy the obvious recurrence relation n= 1 + nk−1. We can now rewrite the definition for the sequence of values maximumk in the form of a definition for the pairs nk, maximumk:

for k = 1 n1 = 1
maximum1 = a[1]

for k > 1 nk = 1 + nk-1
maximumk = MAX(maximumk-1, a[1 + nk -1])

and now the recurrence relations are of the form si = S(si-1), the only —trivial— difference being that in Theorem 8 we started with i = 0 and here with = 1. The trick we called a method shows that we need a second (integer) variable; call it "m". Our state si will associate (with k = i + 1)

maxi = maximumk
mi = n    .

The piece of program now becomes:

max := a[1]; m := 1;
while m < N do begin m := m + 1;
                                   max := MAX(max, a[m])
                          end     .

Again, the order of the two assignment statements is essential.

We have given the above piece of reasoning and the explicit references to the recurrence relation of Theorem 8 because it shows a mechanism leading to the conclusion that the part of the state space on which the repetition operates needs to comprise an additional variable. Even a moderately trained programmer draws this conclusion "intuitively" and from now onwards I shall assume my reader equipped with that intuition. Then —and only then— there is a much shorter path of reasoning that leads to the program we found. It does not consider —"statically" so to speak— the sequence of values s0, s1, … in the sense that it bothers about the values of the subscript i in si. It appeals directly to Theorems 5 and 6 and works in terms of assertions valid (before and after) any execution of the repeatable statement. The price paid for this is the duty to prove termination separately.

Given the base

k = 1 maximum1 = a[1]

and the step

1 < kN maximumk = MAX(maximumk -1, a[k])

the programmer "intuitively" introduces two variables which he calls "maximum" and "k" for short and the relation to be kept invariant is

P : 1 ≤ kN   and   maximum = maximumk     .

(Here the use of "maximum" and "k" stands for the current value of the variable thus named, while "maximumk" stands for the value as given by the recurrence relation. This double use of the same names is tricky but programmers do it. I too.)

The program then consists of two parts: establishing the relation P in accordance with the base and repeatedly increasing k under invariance of relation P, i.e. in accordance with the step.

The initialization

"maximum := a[1]; k := 1"

establishes P (with k = 1), the repetition

while k < N do
begin k := k + 1;
          maximum := MAX(maximum, a[k])

causes the repeatable statement to be executed under the combined relation "B and P ", i.e.

k < N   and   1 ≤ kN   and   maximum = maximumk

which reduces to

1 ≤ k < N   and   maximum = maximumk     .

In order to show that the execution of the repeatable statement under the initial condition (9) leaves relation P valid, it is desirable to distinguish between the values before and after its execution; now it would be confusing to do so with subscripts (why?), therefore we distinguish the values after execution by primes.

Initially we have relation (9); after the assignment k := k + 1 we have the relation k' = k + 1 and from the first part of (9), i.e. 1 ≤ kN, follows 2 ≤ k' ≤ N, which implies

1 ≤ k'    .

The second assignment now becomes effectively maximum := MAX(maximumk, a[k']), resulting in the relation

maximum' = maximumk'     .

Relations (10) and (11) combine to a replica of P, but now for the primed quantities.

Termination follows from the fact that each execution of the repeatable statement involves an effective increase of the integer valued variable k. After termination we have, according to Theorem 5, "P and non B ", i.e.

1 ≤ kN   and   maximum = maximumk    and   non k < N ;

from the first and the last term we conclude k = N and then from the middle part

maximum = maximumN

which concludes the proof.

Exercise. Make a program effectively assigning "prod := X * Y" with integer X and Y, satisfying X ≥ 0, Y ≥ 0
a) using only addition and subtraction
b) using in addition the boolean function "odd(x)", doubling and halving of a number. (The so-called Egyptian multiplication.)

Exercise. Make a program effectively assigning "rem :=REM(X, Y)" with integer X and Y, X ≥ 0, Y > 0, where the function REM(X, Y) is the remainder after the division of X by Y
a) using only addition and subtraction
b) using in addition doubling and halving of a number. Modify both programs in such a way that in addition "quot :=QUOT(X, Y)" will take place. (The so-called Chinese division.)

We conclude this section by an example of the (standard) circumstance in which a recurrence relation should not be translated blindly into a loop. Given two sequences of values

x[1], x[2], … , x[N]        and
y[1], y[2], … , y[N]                    with N ≥ 0     ;

make a program assigning to the boolean variable named "eq" the value true if x[i] = y[i] for all i satisfying 1 ≤ iN and the value false if in that range a value for i exists such that x[i] ≠ y[i]. (The sequence may be empty, in that case "eq" should get the value true.)

How do we define equality for sequences of length N for general N? Again by means of a recurrence relation. Let eqi mean "no difference occurs among the first i pairs"; the sequence of values eqi is given by

for i = 0 eq0 = true
for i > 0 eqi = eqi-1 and x[i] = y[i]     .

The net effect of the program to be made should be eq := eqN.

A blind translation into initialization followed by repetition would lead to

eq := true; i := 0;
while i < N do begin i := i + 1; eq := (eq and x[i] = y[i]) end

Although the above program is correct, the following program, besides being equally correct, is on the average more efficient:

eq := true; i := 0;
while i < N and eq do begin i := i + 1; eq := (x[i] = y[i]) end     .

because it terminates the repetition as soon as a difference has been found.

Exercise. Prove the correctness of the second program.

Back to top

Next chapter: 6. A first example of step-wise program composition