Brun's algorithm

Brun's algorithm is an extension of continued fractions to more than one number. It can be used to find EDOs that contain good approximations to given just ratios.

Algorithm

Say we would like to find EDOs which have good approximations of 5/4, 3/2, and 2/1. First, write the corresponding cents values, largest to smallest, in a vector [math]\displaystyle{ A_{0} }[/math],

\begin{displaymjx} A_{0} = \begin{bmatrix} 1200.00 \\ 701.96 \\ 386.31 \end{bmatrix} \end{displaymjx}

Now write [math]\displaystyle{ XA = A_{0} }[/math], where to begin with [math]\displaystyle{ X }[/math] is the identity matrix and [math]\displaystyle{ A = A_{0} }[/math],

\begin{displaymjx} \begin{bmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} 1200.00 \\ 701.96 \\ 386.31 \end{bmatrix} \\ = \begin{bmatrix} 1200.00 \\ 701.96 \\ 386.31 \end{bmatrix} \end{displaymjx}

Each iteration of the algorithm transforms [math]\displaystyle{ X }[/math] and [math]\displaystyle{ A }[/math] such that

• We preserve [math]\displaystyle{ X A = A_{0} }[/math]
• The elements of [math]\displaystyle{ A }[/math] get smaller
• The first column of [math]\displaystyle{ X }[/math] gives a better EDO approximation of [math]\displaystyle{ A_{0} }[/math]

To do one iteration of the algorithm, firstly subtract the second element from the first element in the cents vector [math]\displaystyle{ A }[/math] and add the first column to the second column in the matrix [math]\displaystyle{ X }[/math], giving

\begin{displaymjx} \begin{bmatrix} 1 & 0 + 1 & 0 \\ 0 & 1 + 0 & 0 \\ 0 & 0 + 0 & 1 \end{bmatrix} \begin{bmatrix} 1200.00 - 701.96 \\ 701.96 \\ 386.31 \end{bmatrix} = \begin{bmatrix} 1 & 1 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} 498.04 \\ 701.96 \\ 386.31 \end{bmatrix} \end{displaymjx}

Secondly, sort the cents vector [math]\displaystyle{ A }[/math] and rearrange the columns of [math]\displaystyle{ X }[/math] in the corresponding order, giving

\begin{displaymjx} \begin{bmatrix} 1 & 1 & 0 \\ 1 & 0 & 0 \\ 0 & 0 & 1 \end{bmatrix} \begin{bmatrix} 701.96 \\ 498.04 \\ 386.31 \end{bmatrix} \end{displaymjx}

Repeating this process gives:

[math]\displaystyle{ X }[/math]	[math]\displaystyle{ A }[/math]	Ratios
[math]\displaystyle{ \begin{bmatrix}1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}1200.00 \\ 701.96 \\ 386.31\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}2/1 \\ 3/2 \\ 5/4\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}1 & 1 & 0 \\ 1 & 0 & 0 \\ 0 & 0 & 1\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}701.96 \\ 498.04 \\ 386.31\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}3/2 \\ 4/3 \\ 5/4\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}2 & 0 & 1 \\ 1 & 0 & 1 \\ 0 & 1 & 0\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}498.04 \\ 386.31 \\ 203.91\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}4/3 \\ 5/4 \\ 9/8\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}2 & 1 & 2 \\ 1 & 1 & 1 \\ 1 & 0 & 0\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}386.31 \\ 203.91 \\ 111.73\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}5/4 \\ 9/8 \\ 16/15\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}3 & 2 & 2 \\ 2 & 1 & 1 \\ 1 & 1 & 0\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}203.91 \\ 182.40 \\ 111.73\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}9/8 \\ 10/9 \\ 16/15\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}5 & 2 & 3 \\ 3 & 1 & 2 \\ 2 & 0 & 1\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}182.40 \\ 111.73 \\ 21.51\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}10/9 \\ 16/15 \\ 81/80\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}7 & 5 & 3 \\ 4 & 3 & 2 \\ 2 & 2 & 1\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}111.73 \\ 70.67 \\ 21.51\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}16/15 \\ 25/24 \\ 81/80\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}12 & 7 & 3 \\ 7 & 4 & 2 \\ 4 & 2 & 1\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}70.67 \\ 41.06 \\ 21.51\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}25/24 \\ 128/125 \\ 81/80\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}19 & 12 & 3 \\ 11 & 7 & 2 \\ 6 & 4 & 1\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}41.06 \\ 29.61 \\ 21.51\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}128/125 \\ 3125/3072 \\ 81/80\end{bmatrix} }[/math]
[math]\displaystyle{ \begin{bmatrix}31 & 3 & 19 \\ 18 & 2 & 11 \\ 10 & 1 & 6\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}29.61 \\ 21.51 \\ 11.45\end{bmatrix} }[/math]	[math]\displaystyle{ \begin{bmatrix}3125/3072 \\ 81/80 \\ 393216/390625\end{bmatrix} }[/math]

Here each row shows one iteration of the algorithm. The Ratios column shows the ratios corresponding to the cents values in [math]\displaystyle{ A }[/math] (from dividing ratios each time we subtract cents values).

The first column of each matrix [math]\displaystyle{ X }[/math] gives an EDO approximation to [math]\displaystyle{ A_{0} }[/math]. For example, the matrix with first column [math]\displaystyle{ [31, 18, 10] }[/math] gives

\begin{displaymjx} \begin{aligned} \frac{1200}{31} \begin{bmatrix} 31 \\ 18 \\ 10 \end{bmatrix} &= \begin{bmatrix} 1200.00 \\ 696.77 \\ 387.10 \end{bmatrix} \\ &\approx \begin{bmatrix} 1200.00 \\ 701.96 \\ 386.31 \end{bmatrix} \end{aligned} \end{displaymjx}

Code

The following Python code implements Brun's algorithm.

import numpy as np

target_ratios = [5/4, 3/2, 2/1]

A = 1200 * np.log2(sorted(target_ratios, reverse=True))
X = np.identity(len(A), dtype=np.int64)

while X.max() <= 100:
    print(f"{X}\n{A}\n")
    A[0] -= A[1]
    X[:, 1] += X[:, 0]
    i = np.argsort(-A)
    A = A[i]
    X = X[:, i]

Other applications

Figure 1: Wilson's scale construction using Brun's algorithm

Erv Wilson used the intervals generated by Brun's algorithm to construct scales. The method starts with a one note scale consisting of the largest interval you are approximating. Then, writing [math]\displaystyle{ a }[/math] and [math]\displaystyle{ b }[/math] for the first two elements of the cents vector [math]\displaystyle{ A }[/math], each time you compute [math]\displaystyle{ a - b }[/math], you split every interval [math]\displaystyle{ a }[/math] in your scale into [math]\displaystyle{ b }[/math] and [math]\displaystyle{ a - b }[/math] (in either order). This gives a sequence of scales in which each contains the previous one. The steps in each scale are given by the elements of the cents vector [math]\displaystyle{ A }[/math] in the corresponding iteration of Brun's algorithm, and their counts are given by the first row of the matrix [math]\displaystyle{ X }[/math]. Figure 1 shows this process for Brun's algorithm applied to 2/1, 3/2, and 5/4.

The Wilsonic iOS app implements scale construction using Brun's method. To choose the splitting order, it uses the heuristic of splitting interval [math]\displaystyle{ a }[/math] into [math]\displaystyle{ [b, a - b] }[/math] unless the previous step was [math]\displaystyle{ b }[/math], in which case it uses the other order [math]\displaystyle{ [a - b, b] }[/math].

Kraig Grady has applied Brun's algorithm analogously to rhythm. He emphasises that the step ordering involves an element of choice:

Unlike Moments of Symmetry patterns, the order in which to place the dividing numbers it is not always clear. Evenness is the key. Best to watch that the shortest duration is flanked by the largest when possible. One ideally wants all composite durations in a row to lie within a close range without overlapping those higher or lower.

Wilson also uses Selmer's algorithm, a variant of Brun's algorithm which subtracts the smallest element of [math]\displaystyle{ A }[/math] from the largest. Wilson calls this the [math]\displaystyle{ a - c }[/math] version of Brun's algorithm since he approximates three ratios and writes [math]\displaystyle{ A = [a, b, c] }[/math]. Selmer's algorithm converges more slowly than Brun's algorithm and so can find some scales that Brun's algorithm does not. For example, starting with the ratios 2/1, 3/2, and 5/4, Brun's algorithm finds scales of sizes 1, 2, 3, 5, 7, 10, 15, 22, 34, 53, ... , while Selmer's algorithm finds sizes 1, 2, 3, 3, 4, 6, 7, 9, 10, 12, 19, 22, 31, 34, 53, ... . Wilson also switches between using [math]\displaystyle{ a - b }[/math] and [math]\displaystyle{ a - c }[/math], a process he calls re-seeding.

Further reading

• Nicolas Chevallier, Best simultaneous Diophantine approximations and multidimensional continued fraction expansions. 2013.

External links

• bruntosaur - Online Brun's algorithm calculator