Binary GCD algorithm

The binary GCD algorithm, also known as Stein's algorithm or the binary Euclidean algorithm, is an algorithm that computes the greatest common divisor of two nonnegative integers. Stein's algorithm uses simpler arithmetic operations than the conventional Euclidean algorithm; it replaces division with arithmetic shifts, comparisons, and subtraction.
Although the algorithm in its contemporary form was first published by the physicist and programmer Josef Stein in 1967, it was known by the 2nd century BCE, in ancient China.

Algorithm

The algorithm finds the GCD of two nonnegative numbers and by repeatedly applying these identities:

: everything divides zero, and is the largest number that divides.
: is a common divisor.
if is odd: is then not a common divisor.
if odd and.

As GCD is commutative, those identities still apply if the operands are swapped:, if is odd, etc.

Implementation

While the above description of the algorithm is mathematically correct, the binary GCD algorithm performs more loop iterations that the Euclidean one, so it only offers a performance advantage if those iterations are substantially faster. Performant software implementations typically differ from it in a few notable ways:

eschewing repeated trial division by in favour of the count trailing zeros primitive and a single bit shift, which is functionally equivalent to repeatedly applying identity 3, but much faster;
expressing the algorithm iteratively rather than recursively: the resulting implementation can be laid out to avoid repeated work, invoking identity 2 at the start and maintaining as invariant that both numbers are odd upon entering the loop, which only needs to implement identities 3 and 4; and
making the loop's body branch-free except for its exit condition : unpredictable branches can have a large, negative impact on performance.

The following is an implementation of the algorithm in C. After the shift at the top of the main loop, the variables and hold the values and. This elimination of the least-significant bit leaves room for the signed difference to be computed on the next line without overflow.

// Count trailing zeros. This is a pedagogical example only;
// for efficiency use __builtin_ctzll or similar.
static int ctz
uint64_t gcd

The following is an implementation of the algorithm in Rust exemplifying those differences, adapted from . The conditional exchange of and compiles down to conditional moves;:

use std::cmp::min;
use std::mem::swap;
pub fn gcd -> u64

Note: The implementation above accepts unsigned integers; given that, the signed case can be handled as follows:

/// Computes the GCD of two signed 64-bit integers
/// The result is unsigned and not always representable as i64: gcd 1 << 63
pub fn signed_gcd -> u64

Complexity

Asymptotically, the algorithm requires steps, where is the number of bits in the larger of the two numbers, as every two steps reduce at least one of the operands by at least a factor of. Each step involves only a few arithmetic operations ; when working with word-sized numbers, each arithmetic operation translates to a single machine operation, so the number of machine operations is on the order of, i.e..
For arbitrarily large numbers, the asymptotic complexity of this algorithm is, as each arithmetic operation involves a linear number of machine operations.
If the numbers can be represented in the machine's memory, i.e. each number's size can be represented by a single machine word, this bound is reduced to:
This is the same as for the Euclidean algorithm, though a more precise analysis by Akhavi and Vallée proved that binary GCD uses about 60% fewer bit operations.

Extensions

The binary GCD algorithm can be extended in several ways, either to output additional information, deal with arbitrarily large integers more efficiently, or to compute GCDs in domains other than the integers.
The extended binary GCD algorithm, analogous to the extended Euclidean algorithm, fits in the first kind of extension, as it provides the Bézout coefficients in addition to the GCD: integers and such that.
In the case of large integers, the best asymptotic complexity is, with the cost of -bit multiplication; this is near-linear and vastly smaller than the binary GCD algorithm's, though concrete implementations only outperform older algorithms for numbers larger than about 64 kilobits. This is achieved by extending the binary GCD algorithm using ideas from the Schönhage–Strassen algorithm for fast integer multiplication.
The binary GCD algorithm has also been extended to domains other than natural numbers, such as Gaussian integers, Eisenstein integers, quadratic rings, and integer rings of number fields.

Historical description

An algorithm for computing the GCD of two numbers was known in ancient China, under the Han dynasty, as a method to reduce fractions:
The phrase "if possible halve it" is ambiguous,

if this applies when either of the numbers become even, the algorithm is the binary GCD algorithm;
if this only applies when both numbers are even, the algorithm is similar to the Euclidean algorithm.