Extended Euclidean Algorithm

Purpose: Alongside $\gcd (a,b)$, find integers $x,y$ with $ax+by=\gcd (a,b)$ (Bézout coefficients), in $O(\log \min (a,b))$.

Algorithm

1. Base case: $\gcd (a,0)=a$, with $x=1$, $y=0$.
2. Otherwise recurse on $(b,a\mathrm{mod}b)$, obtaining $x_1,y_1$ with $b{x}_1+(a\mathrm{mod}b)y_1=g$.
3. Substitute $a\mathrm{mod}b=a-\lfloor a/b\rfloor b$:
   $b{x}_1+\left(a-\lfloor \frac{a}{b}\rfloor b\right)y_1=a{y}_1+b\left(x_1-\lfloor \frac{a}{b}\rfloor y_1\right)=g.$
4. So $x=y_1$ and $y=x_1-\lfloor a/b\rfloor y_1$.

Code

// returns g = gcd(a,b) and sets x, y with a*x + b*y = g
long long extgcd(long long a, long long b, long long &x, long long &y) {
    if (b == 0) { x = 1; y = 0; return a; }
    long long x1, y1;
    long long g = extgcd(b, a % b, x1, y1);
    x = y1;
    y = x1 - (a / b) * y1;
    return g;
}
 
// modular inverse of a mod m (exists iff gcd(a,m) == 1)
long long modinv(long long a, long long m) {
    long long x, y;
    if (extgcd(a, m, x, y) != 1) return -1;   // no inverse
    return ((x % m) + m) % m;
}

An iterative version avoids recursion depth and is what you want inside tight loops:

long long extgcd_it(long long a, long long b, long long &x, long long &y) {
    x = 1; y = 0;
    long long x1 = 0, y1 = 1, a1 = a, b1 = b;
    while (b1) {
        long long q = a1 / b1;
        tie(x, x1) = make_pair(x1, x - q * x1);
        tie(y, y1) = make_pair(y1, y - q * y1);
        tie(a1, b1) = make_pair(b1, a1 - q * b1);
    }
    return a1;
}

Paradigm

Divide and conquer / recursion with back-substitution. The forward pass is the plain Euclidean algorithm; the extension is a linear back-substitution that carries the coefficients up the recursion.

Complexity

• Time: $O(\log \min (a,b))$ — the same Fibonacci-worst-case bound as plain Euclid
• Space: $O(\log \min (a,b))$ recursive, $O(1)$ iterative

Proof of Correctness

By induction on the recursion depth. Base: $a\cdot 1+0\cdot 0=a=\gcd (a,0)$. Step: assume the recursive call returns $x_1,y_1$ with $b{x}_1+(a\mathrm{mod}b)y_1=g$. Substituting $a\mathrm{mod}b=a-\lfloor a/b\rfloor b$ and regrouping (step 3 above) yields exactly $a{y}_1+b(x_1-\lfloor a/b\rfloor y_1)=g$, which is what the algorithm returns. ∎

Size bound: the returned coefficients satisfy $|x|\le b/(2g)$ and $|y|\le a/(2g)$, so they never overflow if $a,b$ fit in the type.

Variants / Use Cases

• Modular inverse — $a^{-1}\mathrm{mod}m$ when $\gcd (a,m)=1$; works for non-prime $m$, unlike Fermat
• Linear Diophantine equations — $ax+by=c$ is solvable iff $g\mid c$; scale the Bézout pair by $c/g$
• CRT for non-coprime moduli — merge congruences using the Bézout coefficients
• Linear congruence $ax\equiv b(\mathrm{mod}m)$
• Continued fractions — the quotients $\lfloor a/b\rfloor$ are exactly the continued fraction expansion of $a/b$
• Garner’s algorithm — mixed-radix CRT reconstruction