Euclidean algorithm

The Euclidean algorithm proceeds in a series of steps, with the output of each step used as the input for the next. Track the steps using an integer counter k, so the initial step corresponds to k = 0, the next step to k = 1, and so on.

Each step begins with two nonnegative remainders r_k−2 and r_k−1, with r_k−2 > r_k−1. The kth step performs division-with-remainder to find the quotient q_k and remainder r_k so that:

${\displaystyle r_{k-2}=q_{k}r_{k-1}+r_{k}\ {\text{ with }}\ r_{k-1}>r_{k}\geq 0.}$

That is, multiples of the smaller number r_k−1 are subtracted from the larger number r_k−2 until the remainder r_k is smaller than r_k−1. Then the algorithm proceeds to the (k+1)th step starting with r_k−1 and r_k.

In the initial step k = 0, the remainders are set to r₋₂ = a and r₋₁ = b, the numbers for which the GCD is sought. In the next step k = 1, the remainders are r₋₁ = b and the remainder r₀ of the initial step, and so on. The algorithm proceeds in a sequence of equations

${\displaystyle {\begin{aligned}a&=q_{0}b+r_{0}\\b&=q_{1}r_{0}+r_{1}\\r_{0}&=q_{2}r_{1}+r_{2}\\r_{1}&=q_{3}r_{2}+r_{3}\\&\,\,\,\vdots \end{aligned}}}$

The algorithm need not be modified if a < b: in that case, the initial quotient is q₀ = 0, the first remainder is r₀ = a, and henceforth r_k−2 > r_k−1 for all k ≥ 1.

Since the remainders are non-negative integers that decrease with every step, the sequence ${\displaystyle r_{-1}>r_{0}>r_{1}>r_{2}>\cdots \geq 0}$ cannot be infinite, so the algorithm must eventually fail to produce the next step; but the division algorithm can always proceed to the (N+1)th step provided r_N > 0. Thus the algorithm must eventually produce a zero remainder r_N = 0.[15] The final nonzero remainder is the greatest common divisor of a and b:

${\displaystyle r_{N-1}=\gcd(a,b).}$

The validity of the Euclidean algorithm can be proven by a two-step argument.[16] In the first step, the final nonzero remainder r_N−1 is shown to divide both a and b. Since it is a common divisor, it must be less than or equal to the greatest common divisor g. In the second step, it is shown that any common divisor of a and b, including g, must divide r_N−1; therefore, g must be less than or equal to r_N−1. These two opposite inequalities imply r_N−1 = g.

To demonstrate that r_N−1 divides both a and b (the first step), r_N−1 divides its predecessor r_N−2

r_N−2 = q_N r_N−1

since the final remainder r_N is zero. r_N−1 also divides its next predecessor r_N−3

r_N−3 = q_N−1 r_N−2 + r_N−1

because it divides both terms on the right-hand side of the equation. Iterating the same argument, r_N−1 divides all the preceding remainders, including a and b. None of the preceding remainders r_N−2, r_N−3, etc. divide a and b, since they leave a remainder. Since r_N−1 is a common divisor of a and b, r_N−1 ≤ g.

In the second step, any natural number c that divides both a and b (in other words, any common divisor of a and b) divides the remainders r_k. By definition, a and b can be written as multiples of c : a = mc and b = nc, where m and n are natural numbers. Therefore, c divides the initial remainder r₀, since r₀ = a − q₀b = mc − q₀nc = (m − q₀n)c. An analogous argument shows that c also divides the subsequent remainders r₁, r₂, etc. Therefore, the greatest common divisor g must divide r_N−1, which implies that g ≤ r_N−1. Since the first part of the argument showed the reverse (r_N−1 ≤ g), it follows that g = r_N−1. Thus, g is the greatest common divisor of all the succeeding pairs:[17][18]

g = gcd(a, b) = gcd(b, r₀) = gcd(r₀, r₁) = … = gcd(r_N−2, r_N−1) = r_N−1.

Subtraction-based animation of the Euclidean algorithm. The initial rectangle has dimensions a = 1071 and b = 462. Squares of size 462×462 are placed within it leaving a 462×147 rectangle. This rectangle is tiled with 147×147 squares until a 21×147 rectangle is left, which in turn is tiled with 21×21 squares, leaving no uncovered area. The smallest square size, 21, is the GCD of 1071 and 462.

For illustration, the Euclidean algorithm can be used to find the greatest common divisor of a = 1071 and b = 462. To begin, multiples of 462 are subtracted from 1071 until the remainder is less than 462. Two such multiples can be subtracted (q₀ = 2), leaving a remainder of 147:

1071 = 2 × 462 + 147.

Then multiples of 147 are subtracted from 462 until the remainder is less than 147. Three multiples can be subtracted (q₁ = 3), leaving a remainder of 21:

462 = 3 × 147 + 21.

Then multiples of 21 are subtracted from 147 until the remainder is less than 21. Seven multiples can be subtracted (q₂ = 7), leaving no remainder:

147 = 7 × 21 + 0.

Since the last remainder is zero, the algorithm ends with 21 as the greatest common divisor of 1071 and 462. This agrees with the gcd(1071, 462) found by prime factorization above. In tabular form, the steps are:

The Euclidean algorithm can be visualized in terms of the tiling analogy given above for the greatest common divisor.[19] Assume that we wish to cover an a×b rectangle with square tiles exactly, where a is the larger of the two numbers. We first attempt to tile the rectangle using b×b square tiles; however, this leaves an r₀×b residual rectangle untiled, where r₀ < b. We then attempt to tile the residual rectangle with r₀×r₀ square tiles. This leaves a second residual rectangle r₁×r₀, which we attempt to tile using r₁×r₁ square tiles, and so on. The sequence ends when there is no residual rectangle, i.e., when the square tiles cover the previous residual rectangle exactly. The length of the sides of the smallest square tile is the GCD of the dimensions of the original rectangle. For example, the smallest square tile in the adjacent figure is 21×21 (shown in red), and 21 is the GCD of 1071 and 462, the dimensions of the original rectangle (shown in green).

At every step k, the Euclidean algorithm computes a quotient q_k and remainder r_k from two numbers r_k−1 and r_k−2

r_k−2 = q_k r_k−1 + r_k

where the r_k is non-negative and is strictly less than the absolute value of r_k−1. The theorem which underlies the definition of the Euclidean division ensures that such a quotient and remainder always exist and are unique.[20]

In Euclid's original version of the algorithm, the quotient and remainder are found by repeated subtraction; that is, r_k−1 is subtracted from r_k−2 repeatedly until the remainder r_k is smaller than r_k−1. After that r_k and r_k−1 are exchanged and the process is iterated. Euclidean division reduces all the steps between two exchanges into a single step, which is thus more efficient. Moreover, the quotients are not needed, thus one may replace Euclidean division by the modulo operation, which gives only the remainder. Thus the iteration of the Euclidean algorithm becomes simply

r_k = r_k−2 mod r_k−1.

Implementations of the algorithm may be expressed in pseudocode. For example, the division-based version may be programmed as[21]

function gcd(a, b)
    while b ≠ 0
        t := b
        b := a mod b
        a := t
    return a

At the beginning of the kth iteration, the variable b holds the latest remainder r_k−1, whereas the variable a holds its predecessor, r_k−2. The step b := a mod b is equivalent to the above recursion formula r_k ≡ r_k−2 mod r_k−1. The temporary variable t holds the value of r_k−1 while the next remainder r_k is being calculated. At the end of the loop iteration, the variable b holds the remainder r_k, whereas the variable a holds its predecessor, r_k−1.

(If negative inputs are allowed, or if the mod function may return negative values, the last line must be changed into return abs(a).)

In the subtraction-based version, which was Euclid's original version, the remainder calculation (b := a mod b) is replaced by repeated subtraction.[22] Contrary to the division-based version, which works with arbitrary integers as input, the subtraction-based version supposes that the input consists of positive integers and stops when a = b:

function gcd(a, b)
    while a ≠ b 
        if a > b
            a := a − b
        else
            b := b − a
    return a

The variables a and b alternate holding the previous remainders r_k−1 and r_k−2. Assume that a is larger than b at the beginning of an iteration; then a equals r_k−2, since r_k−2 > r_k−1. During the loop iteration, a is reduced by multiples of the previous remainder b until a is smaller than b. Then a is the next remainder r_k. Then b is reduced by multiples of a until it is again smaller than a, giving the next remainder r_k+1, and so on.

The recursive version[23] is based on the equality of the GCDs of successive remainders and the stopping condition gcd(r_N−1, 0) = r_N−1.

function gcd(a, b)
    if b = 0
        return a
    else
        return gcd(b, a mod b)

(As above, if negative inputs are allowed, or if the mod function may return negative values, the instruction "return a" must be changed into "return max(a, −a)".)

For illustration, the gcd(1071, 462) is calculated from the equivalent gcd(462, 1071 mod 462) = gcd(462, 147). The latter GCD is calculated from the gcd(147, 462 mod 147) = gcd(147, 21), which in turn is calculated from the gcd(21, 147 mod 21) = gcd(21, 0) = 21.

In another version of Euclid's algorithm, the quotient at each step is increased by one if the resulting negative remainder is smaller in magnitude than the typical positive remainder.[24][25] Previously, the equation

r_k−2 = q_k r_k−1 + r_k

assumed that |r_k−1| > r_k > 0. However, an alternative negative remainder e_k can be computed:

r_k−2 = (q_k + 1) r_k−1 + e_k

if r_k−1 > 0 or

r_k−2 = (q_k − 1) r_k−1 + e_k

if r_k−1 < 0.

If r_k is replaced by e_k. when |e_k| < |r_k|, then one gets a variant of Euclidean algorithm such that

|r_k| ≤ |r_k−1| / 2

at each step.

Leopold Kronecker has shown that this version requires the fewest steps of any version of Euclid's algorithm.[24][25] More generally, it has been proven that, for every input numbers a and b, the number of steps is minimal if and only if q_k is chosen in order that ${\displaystyle \left|{\frac {r_{k+1}}{r_{k}}}\right|<{\frac {1}{\varphi }}\sim 0.618,}$ where ${\displaystyle \varphi }$ is the golden ratio.[26]

Bézout's identity states that the greatest common divisor g of two integers a and b can be represented as a linear sum of the original two numbers a and b.[55] In other words, it is always possible to find integers s and t such that g = sa + tb.[56][57]

The integers s and t can be calculated from the quotients q₀, q₁, etc. by reversing the order of equations in Euclid's algorithm.[58] Beginning with the next-to-last equation, g can be expressed in terms of the quotient q_N−1 and the two preceding remainders, r_N−2 and r_N−3:

g = r_N−1 = r_N−3 − q_N−1 r_N−2 .

Those two remainders can be likewise expressed in terms of their quotients and preceding remainders,

r_N−2 = r_N−4 − q_N−2 r_N−3 and

r_N−3 = r_N−5 − q_N−3 r_N−4 .

Substituting these formulae for r_N−2 and r_N−3 into the first equation yields g as a linear sum of the remainders r_N−4 and r_N−5. The process of substituting remainders by formulae involving their predecessors can be continued until the original numbers a and b are reached:

r₂ = r₀ − q₂ r₁

r₁ = b − q₁ r₀

r₀ = a − q₀ b.

After all the remainders r₀, r₁, etc. have been substituted, the final equation expresses g as a linear sum of a and b, so that g = sa + tb.

The Euclidean algorithm, and thus Bezout's identity, can be generalized to the context of Euclidean domains.

Bézout's identity provides yet another definition of the greatest common divisor g of two numbers a and b.[12] Consider the set of all numbers ua + vb, where u and v are any two integers. Since a and b are both divisible by g, every number in the set is divisible by g. In other words, every number of the set is an integer multiple of g. This is true for every common divisor of a and b. However, unlike other common divisors, the greatest common divisor is a member of the set; by Bézout's identity, choosing u = s and v = t gives g. A smaller common divisor cannot be a member of the set, since every member of the set must be divisible by g. Conversely, any multiple m of g can be obtained by choosing u = ms and v = mt, where s and t are the integers of Bézout's identity. This may be seen by multiplying Bézout's identity by m,

mg = msa + mtb.

Therefore, the set of all numbers ua + vb is equivalent to the set of multiples m of g. In other words, the set of all possible sums of integer multiples of two numbers (a and b) is equivalent to the set of multiples of gcd(a, b). The GCD is said to be the generator of the ideal of a and b. This GCD definition led to the modern abstract algebraic concepts of a principal ideal (an ideal generated by a single element) and a principal ideal domain (a domain in which every ideal is a principal ideal).

Certain problems can be solved using this result.[59] For example, consider two measuring cups of volume a and b. By adding/subtracting u multiples of the first cup and v multiples of the second cup, any volume ua + vb can be measured out. These volumes are all multiples of g = gcd(a, b).

The integers s and t of Bézout's identity can be computed efficiently using the extended Euclidean algorithm. This extension adds two recursive equations to Euclid's algorithm[60]

s_k = s_k−2 − q_ks_k−1

t_k = t_k−2 − q_kt_k−1

with the starting values

s₋₂ = 1, t₋₂ = 0

s₋₁ = 0, t₋₁ = 1.

Using this recursion, Bézout's integers s and t are given by s = s_N and t = t_N, where N+1 is the step on which the algorithm terminates with r_N+1 = 0.

The validity of this approach can be shown by induction. Assume that the recursion formula is correct up to step k − 1 of the algorithm; in other words, assume that

r_j = s_j a + t_j b

for all j less than k. The kth step of the algorithm gives the equation

r_k = r_k−2 − q_kr_k−1.

Since the recursion formula has been assumed to be correct for r_k−2 and r_k−1, they may be expressed in terms of the corresponding s and t variables

r_k = (s_k−2 a + t_k−2 b) − q_k(s_k−1 a + t_k−1 b).

Rearranging this equation yields the recursion formula for step k, as required

r_k = s_k a + t_k b = (s_k−2 − q_ks_k−1) a + (t_k−2 − q_kt_k−1) b.

The integers s and t can also be found using an equivalent matrix method.[61] The sequence of equations of Euclid's algorithm

${\displaystyle {\begin{aligned}a&=q_{0}b+r_{0}\\b&=q_{1}r_{0}+r_{1}\\&\,\,\,\vdots \\r_{N-2}&=q_{N}r_{N-1}+0\end{aligned}}}$

can be written as a product of 2×2 quotient matrices multiplying a two-dimensional remainder vector

${\displaystyle {\begin{pmatrix}a\\b\end{pmatrix}}={\begin{pmatrix}q_{0}&1\\1&0\end{pmatrix}}{\begin{pmatrix}b\\r_{0}\end{pmatrix}}={\begin{pmatrix}q_{0}&1\\1&0\end{pmatrix}}{\begin{pmatrix}q_{1}&1\\1&0\end{pmatrix}}{\begin{pmatrix}r_{0}\\r_{1}\end{pmatrix}}=\cdots =\prod _{i=0}^{N}{\begin{pmatrix}q_{i}&1\\1&0\end{pmatrix}}{\begin{pmatrix}r_{N-1}\\0\end{pmatrix}}\,.}$

Let M represent the product of all the quotient matrices

${\displaystyle \mathbf {M} ={\begin{pmatrix}m_{11}&m_{12}\\m_{21}&m_{22}\end{pmatrix}}=\prod _{i=0}^{N}{\begin{pmatrix}q_{i}&1\\1&0\end{pmatrix}}={\begin{pmatrix}q_{0}&1\\1&0\end{pmatrix}}{\begin{pmatrix}q_{1}&1\\1&0\end{pmatrix}}\cdots {\begin{pmatrix}q_{N}&1\\1&0\end{pmatrix}}\,.}$

This simplifies the Euclidean algorithm to the form

${\displaystyle {\begin{pmatrix}a\\b\end{pmatrix}}=\mathbf {M} {\begin{pmatrix}r_{N-1}\\0\end{pmatrix}}=\mathbf {M} {\begin{pmatrix}g\\0\end{pmatrix}}\,.}$

To express g as a linear sum of a and b, both sides of this equation can be multiplied by the inverse of the matrix M.[61][62] The determinant of M equals (−1)^N+1, since it equals the product of the determinants of the quotient matrices, each of which is negative one. Since the determinant of M is never zero, the vector of the final remainders can be solved using the inverse of M

${\displaystyle {\begin{pmatrix}g\\0\end{pmatrix}}=\mathbf {M} ^{-1}{\begin{pmatrix}a\\b\end{pmatrix}}=(-1)^{N+1}{\begin{pmatrix}m_{22}&-m_{12}\\-m_{21}&m_{11}\end{pmatrix}}{\begin{pmatrix}a\\b\end{pmatrix}}\,.}$

Since the top equation gives

g = (−1)^N+1 ( m₂₂ a − m₁₂ b),

the two integers of Bézout's identity are s = (−1)^N+1m₂₂ and t = (−1)^Nm₁₂. The matrix method is as efficient as the equivalent recursion, with two multiplications and two additions per step of the Euclidean algorithm.

Bézout's identity is essential to many applications of Euclid's algorithm, such as demonstrating the unique factorization of numbers into prime factors.[63] To illustrate this, suppose that a number L can be written as a product of two factors u and v, that is, L = uv. If another number w also divides L but is coprime with u, then w must divide v, by the following argument: If the greatest common divisor of u and w is 1, then integers s and t can be found such that

1 = su + tw

by Bézout's identity. Multiplying both sides by v gives the relation:

v = suv + twv = sL + twv

Since w divides both terms on the right-hand side, it must also divide the left-hand side, v. This result is known as Euclid's lemma.[64] Specifically, if a prime number divides L, then it must divide at least one factor of L. Conversely, if a number w is coprime to each of a series of numbers a₁, a₂, ..., a_n, then w is also coprime to their product, a₁ × a₂ × ... × a_n.[64]

Euclid's lemma suffices to prove that every number has a unique factorization into prime numbers.[65] To see this, assume the contrary, that there are two independent factorizations of L into m and n prime factors, respectively

L = p₁p₂…p_m = q₁q₂…q_n .

Since each prime p divides L by assumption, it must also divide one of the q factors; since each q is prime as well, it must be that p = q. Iteratively dividing by the p factors shows that each p has an equal counterpart q; the two prime factorizations are identical except for their order. The unique factorization of numbers into primes has many applications in mathematical proofs, as shown below.

Plot of a linear Diophantine equation, 9x + 12y = 483. The solutions are shown as blue circles.

Diophantine equations are equations in which the solutions are restricted to integers; they are named after the 3rd-century Alexandrian mathematician Diophantus.[66] A typical linear Diophantine equation seeks integers x and y such that[67]

ax + by = c

where a, b and c are given integers. This can be written as an equation for x in modular arithmetic:

ax ≡ c mod b.

Let g be the greatest common divisor of a and b. Both terms in ax + by are divisible by g; therefore, c must also be divisible by g, or the equation has no solutions. By dividing both sides by c/g, the equation can be reduced to Bezout's identity

sa + tb = g

where s and t can be found by the extended Euclidean algorithm.[68] This provides one solution to the Diophantine equation, x₁ = s (c/g) and y₁ = t (c/g).

In general, a linear Diophantine equation has no solutions, or an infinite number of solutions.[69] To find the latter, consider two solutions, (x₁, y₁) and (x₂, y₂), where

ax₁ + by₁ = c = ax₂ + by₂

or equivalently

a(x₁ − x₂) = b(y₂ − y₁).

Therefore, the smallest difference between two x solutions is b/g, whereas the smallest difference between two y solutions is a/g. Thus, the solutions may be expressed as

x = x₁ − bu/g

y = y₁ + au/g.

By allowing u to vary over all possible integers, an infinite family of solutions can be generated from a single solution (x₁, y₁). If the solutions are required to be positive integers (x > 0, y > 0), only a finite number of solutions may be possible. This restriction on the acceptable solutions allows some systems of Diophantine equations with more unknowns than equations to have a finite number of solutions;[70] this is impossible for a system of linear equations when the solutions can be any real number (see Underdetermined system).

A finite field is a set of numbers with four generalized operations. The operations are called addition, subtraction, multiplication and division and have their usual properties, such as commutativity, associativity and distributivity. An example of a finite field is the set of 13 numbers {0, 1, 2, ..., 12} using modular arithmetic. In this field, the results of any mathematical operation (addition, subtraction, multiplication, or division) is reduced modulo 13; that is, multiples of 13 are added or subtracted until the result is brought within the range 0–12. For example, the result of 5 × 7 = 35 mod 13 = 9. Such finite fields can be defined for any prime p; using more sophisticated definitions, they can also be defined for any power m of a prime p ^m. Finite fields are often called Galois fields, and are abbreviated as GF(p) or GF(p ^m).

In such a field with m numbers, every nonzero element a has a unique modular multiplicative inverse, a⁻¹ such that aa⁻¹ = a⁻¹a ≡ 1 mod m. This inverse can be found by solving the congruence equation ax ≡ 1 mod m,[71] or the equivalent linear Diophantine equation[72]

ax + my = 1.

This equation can be solved by the Euclidean algorithm, as described above. Finding multiplicative inverses is an essential step in the RSA algorithm, which is widely used in electronic commerce; specifically, the equation determines the integer used to decrypt the message.[73] Although the RSA algorithm uses rings rather than fields, the Euclidean algorithm can still be used to find a multiplicative inverse where one exists. The Euclidean algorithm also has other applications in error-correcting codes; for example, it can be used as an alternative to the Berlekamp–Massey algorithm for decoding BCH and Reed–Solomon codes, which are based on Galois fields.[74]

Euclid's algorithm can also be used to solve multiple linear Diophantine equations.[75] Such equations arise in the Chinese remainder theorem, which describes a novel method to represent an integer x. Instead of representing an integer by its digits, it may be represented by its remainders x_i modulo a set of N coprime numbers m_i:[76]

${\displaystyle {\begin{aligned}x_{1}&\equiv x{\pmod {m_{1}}}\\x_{2}&\equiv x{\pmod {m_{2}}}\\&\,\,\,\vdots \\x_{N}&\equiv x{\pmod {m_{N}}}\,.\end{aligned}}}$

The goal is to determine x from its N remainders x_i. The solution is to combine the multiple equations into a single linear Diophantine equation with a much larger modulus M that is the product of all the individual moduli m_i, and define M_i as

${\displaystyle M_{i}={\frac {M}{m_{i}}}.}$

Thus, each M_i is the product of all the moduli except m_i. The solution depends on finding N new numbers h_i such that

${\displaystyle M_{i}h_{i}\equiv 1{\pmod {m_{i}}}\,.}$

With these numbers h_i, any integer x can be reconstructed from its remainders x_i by the equation

${\displaystyle x\equiv (x_{1}M_{1}h_{1}+x_{2}M_{2}h_{2}+\cdots +x_{N}M_{N}h_{N}){\pmod {M}}\,.}$

Since these numbers h_i are the multiplicative inverses of the M_i, they may be found using Euclid's algorithm as described in the previous subsection.

The Euclidean algorithm can be used to arrange the set of all positive rational numbers into an infinite binary search tree, called the Stern–Brocot tree. The number 1 (expressed as a fraction 1/1) is placed at the root of the tree, and the location of any other number a/b can be found by computing gcd(a,b) using the original form of the Euclidean algorithm, in which each step replaces the larger of the two given numbers by its difference with the smaller number (not its remainder), stopping when two equal numbers are reached. A step of the Euclidean algorithm that replaces the first of the two numbers corresponds to a step in the tree from a node to its right child, and a step that replaces the second of the two numbers corresponds to a step in the tree from a node to its left child. The sequence of steps constructed in this way does not depend on whether a/b is given in lowest terms, and forms a path from the root to a node containing the number a/b.[77] This fact can be used to prove that each positive rational number appears exactly once in this tree.

For example, 3/4 can be found by starting at the root, going to the left once, then to the right twice:

The Stern–Brocot tree, and the Stern–Brocot sequences of order i for i = 1, 2, 3, 4

${\displaystyle {\begin{aligned}&\gcd(3,4)&\leftarrow \\={}&\gcd(3,1)&\rightarrow \\={}&\gcd(2,1)&\rightarrow \\={}&\gcd(1,1).\end{aligned}}}$

The Euclidean algorithm has almost the same relationship to another binary tree on the rational numbers called the Calkin–Wilf tree. The difference is that the path is reversed: instead of producing a path from the root of the tree to a target, it produces a path from the target to the root.

The Euclidean algorithm has a close relationship with continued fractions.[78] The sequence of equations can be written in the form

${\displaystyle {\begin{aligned}{\frac {a}{b}}&=q_{0}+{\frac {r_{0}}{b}}\\{\frac {b}{r_{0}}}&=q_{1}+{\frac {r_{1}}{r_{0}}}\\{\frac {r_{0}}{r_{1}}}&=q_{2}+{\frac {r_{2}}{r_{1}}}\\&\,\,\,\vdots \\{\frac {r_{k-2}}{r_{k-1}}}&=q_{k}+{\frac {r_{k}}{r_{k-1}}}\\&\,\,\,\vdots \\{\frac {r_{N-2}}{r_{N-1}}}&=q_{N}\,.\end{aligned}}}$

The last term on the right-hand side always equals the inverse of the left-hand side of the next equation. Thus, the first two equations may be combined to form

${\displaystyle {\frac {a}{b}}=q_{0}+{\cfrac {1}{q_{1}+{\cfrac {r_{1}}{r_{0}}}}}\,.}$

The third equation may be used to substitute the denominator term r₁/r₀, yielding

${\displaystyle {\frac {a}{b}}=q_{0}+{\cfrac {1}{q_{1}+{\cfrac {1}{q_{2}+{\cfrac {r_{2}}{r_{1}}}}}}}\,.}$

The final ratio of remainders r_k/r_k−1 can always be replaced using the next equation in the series, up to the final equation. The result is a continued fraction

${\displaystyle {\frac {a}{b}}=q_{0}+{\cfrac {1}{q_{1}+{\cfrac {1}{q_{2}+{\cfrac {1}{\ddots +{\cfrac {1}{q_{N}}}}}}}}}=[q_{0};q_{1},q_{2},\ldots ,q_{N}]\,.}$

In the worked example above, the gcd(1071, 462) was calculated, and the quotients q_k were 2, 3 and 7, respectively. Therefore, the fraction 1071/462 may be written

${\displaystyle {\frac {1071}{462}}=2+{\cfrac {1}{3+{\cfrac {1}{7}}}}=[2;3,7]}$

as can be confirmed by calculation.

Calculating a greatest common divisor is an essential step in several integer factorization algorithms,[79] such as Pollard's rho algorithm,[80] Shor's algorithm,[81] Dixon's factorization method[82] and the Lenstra elliptic curve factorization.[83] The Euclidean algorithm may be used to find this GCD efficiently. Continued fraction factorization uses continued fractions, which are determined using Euclid's algorithm.[84]

The number of steps to calculate the GCD of two natural numbers, a and b, may be denoted by T(a, b).[95] If g is the GCD of a and b, then a = mg and b = ng for two coprime numbers m and n. Then

T(a, b) = T(m, n)

as may be seen by dividing all the steps in the Euclidean algorithm by g.[96] By the same argument, the number of steps remains the same if a and b are multiplied by a common factor w: T(a, b) = T(wa, wb). Therefore, the number of steps T may vary dramatically between neighboring pairs of numbers, such as T(a, b) and T(a, b + 1), depending on the size of the two GCDs.

The recursive nature of the Euclidean algorithm gives another equation

T(a, b) = 1 + T(b, r₀) = 2 + T(r₀, r₁) = … = N + T(r_N−2, r_N−1) = N + 1

where T(x, 0) = 0 by assumption.[95]

If the Euclidean algorithm requires N steps for a pair of natural numbers a > b > 0, the smallest values of a and b for which this is true are the Fibonacci numbers F_N+2 and F_N+1, respectively.[97] More precisely, if the Euclidean algorithm requires N steps for the pair a > b, then one has a ≥ F_N+2 and b ≥ F_N+1. This can be shown by induction.[98] If N = 1, b divides a with no remainder; the smallest natural numbers for which this is true is b = 1 and a = 2, which are F₂ and F₃, respectively. Now assume that the result holds for all values of N up to M − 1. The first step of the M-step algorithm is a = q₀b + r₀, and the Euclidean algorithm requires M − 1 steps for the pair b > r₀. By induction hypothesis, one has b ≥ F_M+1 and r₀ ≥ F_M. Therefore, a = q₀b + r₀ ≥ b + r₀ ≥ F_M+1 + F_M = F_M+2, which is the desired inequality. This proof, published by Gabriel Lamé in 1844, represents the beginning of computational complexity theory,[99] and also the first practical application of the Fibonacci numbers.[97]

This result suffices to show that the number of steps in Euclid's algorithm can never be more than five times the number of its digits (base 10).[100] For if the algorithm requires N steps, then b is greater than or equal to F_N+1 which in turn is greater than or equal to φ^N−1, where φ is the golden ratio. Since b ≥ φ^N−1, then N − 1 ≤ log_φb. Since log₁₀φ > 1/5, (N − 1)/5 < log₁₀φ log_φb = log₁₀b. Thus, N ≤ 5 log₁₀b. Thus, the Euclidean algorithm always needs less than O(h) divisions, where h is the number of digits in the smaller number b.

The average number of steps taken by the Euclidean algorithm has been defined in three different ways. The first definition is the average time T(a) required to calculate the GCD of a given number a and a smaller natural number b chosen with equal probability from the integers 0 to a − 1[95]

${\displaystyle T(a)={\frac {1}{a}}\sum _{0\leq b<a}T(a,b).}$

However, since T(a, b) fluctuates dramatically with the GCD of the two numbers, the averaged function T(a) is likewise "noisy".[101]

To reduce this noise, a second average τ(a) is taken over all numbers coprime with a

${\displaystyle \tau (a)={\frac {1}{\varphi (a)}}\sum _{\begin{smallmatrix}0\leq b<a\\\gcd(a,b)=1\end{smallmatrix}}T(a,b).}$

There are φ(a) coprime integers less than a, where φ is Euler's totient function. This tau average grows smoothly with a[102][103]

${\displaystyle \tau (a)={\frac {12}{\pi ^{2}}}\ln 2\ln a+C+O(a^{-1/6-\varepsilon })}$

with the residual error being of order a^{−(1/6) + ε}, where ε is infinitesimal. The constant C in this formula is called Porter's constant[104] and equals

${\displaystyle C=-{\frac {1}{2}}+{\frac {6\ln 2}{\pi ^{2}}}\left(4\gamma -{\frac {24}{\pi ^{2}}}\zeta '(2)+3\ln 2-2\right)\approx 1.467}$

where γ is the Euler–Mascheroni constant and ζ' is the derivative of the Riemann zeta function.[105][106] The leading coefficient (12/π²) ln 2 was determined by two independent methods.[107][108]

Since the first average can be calculated from the tau average by summing over the divisors d of a[109]

${\displaystyle T(a)={\frac {1}{a}}\sum _{d\mid a}\varphi (d)\tau (d)}$

it can be approximated by the formula[110]

${\displaystyle T(a)\approx C+{\frac {12}{\pi ^{2}}}\ln 2\left(\ln a-\sum _{d\mid a}{\frac {\Lambda (d)}{d}}\right)}$

where Λ(d) is the Mangoldt function.[111]

A third average Y(n) is defined as the mean number of steps required when both a and b are chosen randomly (with uniform distribution) from 1 to n[110]

${\displaystyle Y(n)={\frac {1}{n^{2}}}\sum _{a=1}^{n}\sum _{b=1}^{n}T(a,b)={\frac {1}{n}}\sum _{a=1}^{n}T(a).}$

Substituting the approximate formula for T(a) into this equation yields an estimate for Y(n)[112]

${\displaystyle Y(n)\approx {\frac {12}{\pi ^{2}}}\ln 2\ln n+0.06.}$

In each step k of the Euclidean algorithm, the quotient q_k and remainder r_k are computed for a given pair of integers r_k−2 and r_k−1

r_k−2 = q_k r_k−1 + r_k.

The computational expense per step is associated chiefly with finding q_k, since the remainder r_k can be calculated quickly from r_k−2, r_k−1, and q_k

r_k = r_k−2 − q_k r_k−1.

The computational expense of dividing h-bit numbers scales as O(h(ℓ+1)), where ℓ is the length of the quotient.[92]

For comparison, Euclid's original subtraction-based algorithm can be much slower. A single integer division is equivalent to the quotient q number of subtractions. If the ratio of a and b is very large, the quotient is large and many subtractions will be required. On the other hand, it has been shown that the quotients are very likely to be small integers. The probability of a given quotient q is approximately ln |u/(u − 1)| where u = (q + 1)².[113] For illustration, the probability of a quotient of 1, 2, 3, or 4 is roughly 41.5%, 17.0%, 9.3%, and 5.9%, respectively. Since the operation of subtraction is faster than division, particularly for large numbers,[114] the subtraction-based Euclid's algorithm is competitive with the division-based version.[115] This is exploited in the binary version of Euclid's algorithm.[116]

Combining the estimated number of steps with the estimated computational expense per step shows that the Euclid's algorithm grows quadratically (h²) with the average number of digits h in the initial two numbers a and b. Let h₀, h₁, ..., h_N−1 represent the number of digits in the successive remainders r₀, r₁, ..., r_N−1. Since the number of steps N grows linearly with h, the running time is bounded by

${\displaystyle O{\Big (}\sum _{i<N}h_{i}(h_{i}-h_{i+1}+2){\Big )}\subseteq O{\Big (}h\sum _{i<N}(h_{i}-h_{i+1}+2){\Big )}\subseteq O(h(h_{0}+2N))\subseteq O(h^{2}).}$

Euclid's algorithm is widely used in practice, especially for small numbers, due to its simplicity.[117] For comparison, the efficiency of alternatives to Euclid's algorithm may be determined.

One inefficient approach to finding the GCD of two natural numbers a and b is to calculate all their common divisors; the GCD is then the largest common divisor. The common divisors can be found by dividing both numbers by successive integers from 2 to the smaller number b. The number of steps of this approach grows linearly with b, or exponentially in the number of digits. Another inefficient approach is to find the prime factors of one or both numbers. As noted above, the GCD equals the product of the prime factors shared by the two numbers a and b.[8] Present methods for prime factorization are also inefficient; many modern cryptography systems even rely on that inefficiency.[11]

The binary GCD algorithm is an efficient alternative that substitutes division with faster operations by exploiting the binary representation used by computers.[118][119] However, this alternative also scales like O(h²). It is generally faster than the Euclidean algorithm on real computers, even though it scales in the same way.[93] Additional efficiency can be gleaned by examining only the leading digits of the two numbers a and b.[120][121] The binary algorithm can be extended to other bases (k-ary algorithms),[122] with up to fivefold increases in speed.[123] Lehmer's GCD algorithm uses the same general principle as the binary algorithm to speed up GCD computations in arbitrary bases.

A recursive approach for very large integers (with more than 25,000 digits) leads to quasilinear integer GCD algorithms,[124] such as those of Schönhage,[125][126] and Stehlé and Zimmermann.[127] These algorithms exploit the 2×2 matrix form of the Euclidean algorithm given above. These quasilinear methods generally scale as O(h log h² log log h).[93][94]

Euclid's algorithm can be applied to real numbers, as described by Euclid in Book 10 of his Elements. The goal of the algorithm is to identify a real number g such that two given real numbers, a and b, are integer multiples of it: a = mg and b = ng, where m and n are integers.[28] This identification is equivalent to finding an integer relation among the real numbers a and b; that is, it determines integers s and t such that sa + tb = 0. If such an equation is possible, a and b are called commensurable lengths, otherwise they are incommensurable lengths.[131][132]

The real-number Euclidean algorithm differs from its integer counterpart in two respects. First, the remainders r_k are real numbers, although the quotients q_k are integers as before. Second, the algorithm is not guaranteed to end in a finite number N of steps. If it does, the fraction a/b is a rational number, i.e., the ratio of two integers

${\displaystyle {\frac {a}{b}}={\frac {mg}{ng}}={\frac {m}{n}},}$

and can be written as a finite continued fraction [q₀; q₁, q₂, ..., q_N]. If the algorithm does not stop, the fraction a/b is an irrational number and can be described by an infinite continued fraction [q₀; q₁, q₂, …].[133] Examples of infinite continued fractions are the golden ratio φ = [1; 1, 1, ...] and the square root of two, √2 = [1; 2, 2, ...].[134] The algorithm is unlikely to stop, since almost all ratios a/b of two real numbers are irrational.[135]

An infinite continued fraction may be truncated at a step k [q₀; q₁, q₂, ..., q_k] to yield an approximation to a/b that improves as k is increased. The approximation is described by convergents m_k/n_k; the numerator and denominators are coprime and obey the recurrence relation

${\displaystyle {\begin{aligned}m_{k}&=q_{k}m_{k-1}+m_{k-2}\\n_{k}&=q_{k}n_{k-1}+n_{k-2},\end{aligned}}}$

where m₋₁ = n₋₂ = 1 and m₋₂ = n₋₁ = 0 are the initial values of the recursion. The convergent m_k/n_k is the best rational number approximation to a/b with denominator n_k:[136]

${\displaystyle \left|{\frac {a}{b}}-{\frac {m_{k}}{n_{k}}}\right|<{\frac {1}{n_{k}^{2}}}.}$

Polynomials in a single variable x can be added, multiplied and factored into irreducible polynomials, which are the analogs of the prime numbers for integers. The greatest common divisor polynomial g(x) of two polynomials a(x) and b(x) is defined as the product of their shared irreducible polynomials, which can be identified using the Euclidean algorithm.[128] The basic procedure is similar to that for integers. At each step k, a quotient polynomial q_k(x) and a remainder polynomial r_k(x) are identified to satisfy the recursive equation

${\displaystyle r_{k-2}(x)=q_{k}(x)r_{k-1}(x)+r_{k}(x),}$

where r₋₂(x) = a(x) and r₋₁(x) = b(x). Each quotient polynomial is chosen such that each remainder is either zero or has a degree that is smaller than the degree of its predecessor: deg[r_k(x)] < deg[r_k−1(x)]. Since the degree is a nonnegative integer, and since it decreases with every step, the Euclidean algorithm concludes in a finite number of steps. The last nonzero remainder is the greatest common divisor of the original two polynomials, a(x) and b(x).[137]

For example, consider the following two quartic polynomials, which each factor into two quadratic polynomials

${\displaystyle {\begin{aligned}a(x)&=x^{4}-4x^{3}+4x^{2}-3x+14=(x^{2}-5x+7)(x^{2}+x+2)\qquad {\text{and}}\\b(x)&=x^{4}+8x^{3}+12x^{2}+17x+6=(x^{2}+7x+3)(x^{2}+x+2).\end{aligned}}}$

Dividing a(x) by b(x) yields a remainder r₀(x) = x³ + (2/3)x² + (5/3)x − (2/3). In the next step, b(x) is divided by r₀(x) yielding a remainder r₁(x) = x² + x + 2. Finally, dividing r₀(x) by r₁(x) yields a zero remainder, indicating that r₁(x) is the greatest common divisor polynomial of a(x) and b(x), consistent with their factorization.

Many of the applications described above for integers carry over to polynomials.[138] The Euclidean algorithm can be used to solve linear Diophantine equations and Chinese remainder problems for polynomials; continued fractions of polynomials can also be defined.

The polynomial Euclidean algorithm has other applications, such as Sturm chains, a method for counting the zeros of a polynomial that lie inside a given real interval.[139] This in turn has applications in several areas, such as the Routh–Hurwitz stability criterion in control theory.[140]

Finally, the coefficients of the polynomials need not be drawn from integers, real numbers or even the complex numbers. For example, the coefficients may be drawn from a general field, such as the finite fields GF(p) described above. The corresponding conclusions about the Euclidean algorithm and its applications hold even for such polynomials.[128]

Distribution of Gaussian primes u + vi in the complex plane, with norms u² + v² less than 500

The Gaussian integers are complex numbers of the form α = u + vi, where u and v are ordinary integers[note 2] and i is the square root of negative one.[141] By defining an analog of the Euclidean algorithm, Gaussian integers can be shown to be uniquely factorizable, by the argument above.[42] This unique factorization is helpful in many applications, such as deriving all Pythagorean triples or proving Fermat's theorem on sums of two squares.[141] In general, the Euclidean algorithm is convenient in such applications, but not essential; for example, the theorems can often be proven by other arguments.

The Euclidean algorithm developed for two Gaussian integers α and β is nearly the same as that for ordinary integers,[142] but differs in two respects. As before, we set r₋₂ = α and r₋₁ = β, and the task at each step k is to identify a quotient q_k and a remainder r_k such that

${\displaystyle r_{k}=r_{k-2}-q_{k}r_{k-1},}$

where every remainder is strictly smaller than its predecessor: |r_k| < |r_k−1|. The first difference is that the quotients and remainders are themselves Gaussian integers, and thus are complex numbers. The quotients q_k are generally found by rounding the real and complex parts of the exact ratio (such as the complex number α/β) to the nearest integers.[142] The second difference lies in the necessity of defining how one complex remainder can be "smaller" than another. To do this, a norm function f(u + vi) = u² + v² is defined, which converts every Gaussian integer u + vi into an ordinary integer. After each step k of the Euclidean algorithm, the norm of the remainder f(r_k) is smaller than the norm of the preceding remainder, f(r_k−1). Since the norm is a nonnegative integer and decreases with every step, the Euclidean algorithm for Gaussian integers ends in a finite number of steps.[143] The final nonzero remainder is gcd(α, β), the Gaussian integer of largest norm that divides both α and β; it is unique up to multiplication by a unit, ±1 or ±i.[144]

Many of the other applications of the Euclidean algorithm carry over to Gaussian integers. For example, it can be used to solve linear Diophantine equations and Chinese remainder problems for Gaussian integers;[145] continued fractions of Gaussian integers can also be defined.[142]

A set of elements under two binary operations, denoted as addition and multiplication, is called a Euclidean domain if it forms a commutative ring R and, roughly speaking, if a generalized Euclidean algorithm can be performed on them.[146][147] The two operations of such a ring need not be the addition and multiplication of ordinary arithmetic; rather, they can be more general, such as the operations of a mathematical group or monoid. Nevertheless, these general operations should respect many of the laws governing ordinary arithmetic, such as commutativity, associativity and distributivity.

The generalized Euclidean algorithm requires a Euclidean function, i.e., a mapping f from R into the set of nonnegative integers such that, for any two nonzero elements a and b in R, there exist q and r in R such that a = qb + r and f(r) < f(b).[148] Examples of such mappings are the absolute value for integers, the degree for univariate polynomials, and the norm for Gaussian integers above.[149][150] The basic principle is that each step of the algorithm reduces f inexorably; hence, if f can be reduced only a finite number of times, the algorithm must stop in a finite number of steps. This principle relies on the well-ordering property of the non-negative integers, which asserts that every non-empty set of non-negative integers has a smallest member.[151]

The fundamental theorem of arithmetic applies to any Euclidean domain: Any number from a Euclidean domain can be factored uniquely into irreducible elements. Any Euclidean domain is a unique factorization domain (UFD), although the converse is not true.[151] The Euclidean domains and the UFD's are subclasses of the GCD domains, domains in which a greatest common divisor of two numbers always exists.[152] In other words, a greatest common divisor may exist (for all pairs of elements in a domain), although it may not be possible to find it using a Euclidean algorithm. A Euclidean domain is always a principal ideal domain (PID), an integral domain in which every ideal is a principal ideal.[153] Again, the converse is not true: not every PID is a Euclidean domain.

The unique factorization of Euclidean domains is useful in many applications. For example, the unique factorization of the Gaussian integers is convenient in deriving formulae for all Pythagorean triples and in proving Fermat's theorem on sums of two squares.[141] Unique factorization was also a key element in an attempted proof of Fermat's Last Theorem published in 1847 by Gabriel Lamé, the same mathematician who analyzed the efficiency of Euclid's algorithm, based on a suggestion of Joseph Liouville.[154] Lamé's approach required the unique factorization of numbers of the form x + ωy, where x and y are integers, and ω = e^2iπ/n is an nth root of 1, that is, ωⁿ = 1. Although this approach succeeds for some values of n (such as n = 3, the Eisenstein integers), in general such numbers do not factor uniquely. This failure of unique factorization in some cyclotomic fields led Ernst Kummer to the concept of ideal numbers and, later, Richard Dedekind to ideals.[155]

Distribution of Eisenstein primes u + vω in the complex plane, with norms less than 500. The number ω is a cube root of 1.

The quadratic integer rings are helpful to illustrate Euclidean domains. Quadratic integers are generalizations of the Gaussian integers in which the imaginary unit i is replaced by a number ω. Thus, they have the form u + vω, where u and v are integers and ω has one of two forms, depending on a parameter D. If D does not equal a multiple of four plus one, then

${\displaystyle \omega ={\sqrt {D}}.}$

If, however, D does equal a multiple of four plus one, then

${\displaystyle \omega ={\frac {1+{\sqrt {D}}}{2}}.}$

If the function f corresponds to a norm function, such as that used to order the Gaussian integers above, then the domain is known as norm-Euclidean. The norm-Euclidean rings of quadratic integers are exactly those where D is one of the values −11, −7, −3, −2, −1, 2, 3, 5, 6, 7, 11, 13, 17, 19, 21, 29, 33, 37, 41, 57, or 73.[156][157] The cases D = −1 and D = −3 yield the Gaussian integers and Eisenstein integers, respectively.

If f is allowed to be any Euclidean function, then the list of possible values of D for which the domain is Euclidean is not yet known.[158] The first example of a Euclidean domain that was not norm-Euclidean (with D = 69) was published in 1994.[158] In 1973, Weinberger proved that a quadratic integer ring with D > 0 is Euclidean if, and only if, it is a principal ideal domain, provided that the generalized Riemann hypothesis holds.[129]

The Euclidean algorithm may be applied to some noncommutative rings such as the set of Hurwitz quaternions.[130] Let α and β represent two elements from such a ring. They have a common right divisor δ if α = ξδ and β = ηδ for some choice of ξ and η in the ring. Similarly, they have a common left divisor if α = dξ and β = dη for some choice of ξ and η in the ring. Since multiplication is not commutative, there are two versions of the Euclidean algorithm, one for right divisors and one for left divisors.[130] Choosing the right divisors, the first step in finding the gcd(α, β) by the Euclidean algorithm can be written

${\displaystyle \rho _{0}=\alpha -\psi _{0}\beta =(\xi -\psi _{0}\eta )\delta ,}$

where ψ₀ represents the quotient and ρ₀ the remainder. This equation shows that any common right divisor of α and β is likewise a common divisor of the remainder ρ₀. The analogous equation for the left divisors would be

${\displaystyle \rho _{0}=\alpha -\beta \psi _{0}=\delta (\xi -\eta \psi _{0}).}$

With either choice, the process is repeated as above until the greatest common right or left divisor is identified. As in the Euclidean domain, the "size" of the remainder ρ₀ (formally, its norm) must be strictly smaller than β, and there must be only a finite number of possible sizes for ρ₀, so that the algorithm is guaranteed to terminate.[159]

Most of the results for the GCD carry over to noncommutative numbers. For example, Bézout's identity states that the right gcd(α, β) can be expressed as a linear combination of α and β.[160] In other words, there are numbers σ and τ such that

${\displaystyle \Gamma _{\text{right}}=\sigma \alpha +\tau \beta .}$

The analogous identity for the left GCD is nearly the same:

${\displaystyle \Gamma _{\text{left}}=\alpha \sigma +\beta \tau .}$

Bézout's identity can be used to solve Diophantine equations. For instance, one of the standard proofs of Lagrange's four-square theorem, that every positive integer can be represented as a sum of four squares, is based on quaternion GCDs in this way.[159]