Open main menu

Wikibooks β

Algorithm Implementation/Mathematics/Extended Euclidean algorithm

< Algorithm Implementation‎ | Mathematics

Contents

CEdit

Recursive algorithmEdit

//https://en.wikipedia.org/wiki/Extended_Euclidean_algorithm
//clrs p.937

#include <stdio.h>
#include <stdlib.h>

int mod (int a, int b){
	return a %b; 
}	

int *extendedEuclid (int a, int b){
	int *dxy = (int *)malloc(sizeof(int) *3);

	if (b ==0){
		dxy[0] =a; dxy[1] =1; dxy[2] =0;
		
		return dxy;
	}
	else{
		int t, t2;
		dxy = extendedEuclid(b, mod(a, b));
		t =dxy[1];
		t2 =dxy[2];
		dxy[1] =dxy[2];
		dxy[2] = t - a/b *t2;

		return dxy;
	}
}

int main(void)
{
	int a =99, b =78;
	int *ptr;

	ptr =extendedEuclid (a, b);
	printf("%d = %d * %d + %d * %d \n",ptr[0], a, ptr[1], b, ptr[2] );

	return 0;		
}

PythonEdit

Both functions take positive integers a, b as input, and return a triple (g, x, y), such that ax + by = g = gcd(a, b).

Iterative algorithmEdit

def xgcd(b, n):
    x0, x1, y0, y1 = 1, 0, 0, 1
    while n != 0:
        q, b, n = b // n, n, b % n
        x0, x1 = x1, x0 - q * x1
        y0, y1 = y1, y0 - q * y1
    return  b, x0, y0

Recursive algorithmEdit

import sys
sys.setrecursionlimit(1000000)  # long type,32bit OS 4B,64bit OS 8B(1bit for sign)

# return (g, x, y) a*x + b*y = gcd(x, y)
def egcd(a, b):
    if a == 0:
        return (b, 0, 1)
    else:
        g, x, y = egcd(b % a, a)
        return (g, y - (b // a) * x, x)

Modular inverseEdit

An application of extended GCD algorithm to finding modular inverses:

# x = mulinv(b) mod n, (x * b) % n == 1
def mulinv(b, n):
    g, x, _ = egcd(b, n)
    if g == 1:
        return x % n

HaskellEdit

UnextendedEdit

euclid :: Integral a => a -> a -> a
euclid 0 b = abs b
euclid a 0 = abs a
euclid a b = euclid b $ rem a b

ExtendedEdit

eGCD :: Integer -> Integer -> (Integer,Integer,Integer)
eGCD 0 b = (b, 0, 1)
eGCD a b = let (g, s, t) = eGCD (b `mod` a) a
           in (g, t - (b `div` a) * s, s)