经典算法（4）- 用欧几里得算法实现扩展的最大公约数(Extended GCD)

扩展的gcd算法即除了计算gcd(m,n)还要计算整数x和y，使之满足gcd(m,n) = m.x + n.y。 下面的算法中使用迭代方式。 extendedGCD2方法是extendedGCD的简化版本，考虑到在初值向量r{-1} = [1 0], r{0} = [0 1]下，满足递推关系：r{i} = r{i-2} - q{i}.r{i-1}。

采用Euclid's算法时，不仅要r（余数）的值，还需要q(商)的值。

本例实现参考了Wikipedia中介绍的迭代方法：http://en.wikipedia.org/wiki/Extended_Euclidean_algorithm

/**
 * 
 * @author ljs 2011-5-17
 * 
 * solve extended gcd using Euclid's Algorithm and Iterative Method
 *
 */
public class ExtGCD_Euclid {	
	/**
	 * caculate x,y to satisfy Bezout Identity: m.x + n.y = gcd(m,n)
	 * @param m
	 * @param n
	 * @return {d,x,y} for: d = m.x + n.y;
	 */
	public static int[] extendedGCD(int m,int n){
		m = (m<0)?-m:m;
		n = (n<0)?-n:n;
		
		if(n == 0 ){
			return new int[]{m,1,0};
		}
		
		int r = m % n;
		int q = m / n;
		if(r ==0)
			return new int[]{n,0,1};
		else{
			m=n;
			n=r;
			r = m % n;
			int lastq = q;
			q = m / n;
			if(r==0){
				return new int[]{n,1,-lastq};
			}else{
				int coeffM1=1;
				int coeffN1=-lastq;
				//r2 = n - q2.r1
				int coeffM2=-q;
				int coeffN2=1+lastq * q;
				
				while(n!=0){
					m=n;
					n=r;
					r = m%n;
					q = m / n;
					if(r==0){
						break;
					}else{			
						//for k-th loop, r{k} = r{k-2} - q{k}.r{k-1}, here q{k} is quotient = m{k} / n{k}			
						int coeffM3=coeffM1 - q*coeffM2;
						int coeffN3=coeffN1 - q*coeffN2;	
						coeffM1 = coeffM2;
						coeffN1 = coeffN2;
						coeffM2 = coeffM3;
						coeffN2 = coeffN3;
					}						
				}
				return new int[]{n,coeffM2,coeffN2};
			}
		}		
	}
	
	/**
	 * simplified from the above extendedGCD(m,n)
	 * @see also "http://en.wikipedia.org/wiki/Extended_Euclidean_algorithm"
	 * @param m
	 * @param n
	 * @return {d,x,y} for: d = m.x + n.y;
	 */
	public static int[] extendedGCD2(int m,int n){
		m = (m<0)?-m:m;
		n = (n<0)?-n:n;
		
		int r00=1,r01=0;  //[1 0]
		int r10=0,r11=1;  //[0 1]
		
		int d = 0;
		int x = r00;
		int y = r01;
		while(n!=0){
			int remainder = m%n;
			int quotient = m/n;
			
			if(remainder == 0){			
				x = r10;
				y = r11;
			}else{
				int tmp0 = r00 - quotient*r10;
				int tmp1 = r01 - quotient*r11;
				r00=r10;
				r01=r11;
				r10 = tmp0;
				r11 = tmp1;
			}			
			
			m = n;
			n = remainder;
		}
		d = m;
		
		return new int[]{d,x,y};
	}
	

	public static void print(int m,int n,int[] extGCDResult){
		m = (m<0)?-m:m;
		n = (n<0)?-n:n;
		System.out.format("extended gcd of %d and %d is: %d = %d.{%d} + %d.{%d}%5s%n",m,n,extGCDResult[0],m,extGCDResult[1],n,extGCDResult[2],(m*extGCDResult[1] + n * extGCDResult[2] == extGCDResult[0])?"OK":"WRONG!!!");		
	}
	
	public static void main(String[] args) {
		
		int m = -18;
		int n= 12;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		//co-prime
		m = 15;
		n= 28;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
				
		m = 6;
		n= 3;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 6;
		n= 3;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 6;
		n= 0;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 0;
		n= 6;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 0;
		n= 0;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 1;
		n= 1;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 3;
		n= 3;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 2;
		n= 2;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 1;
		n= 4;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 4;
		n= 1;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 10;
		n= 14;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 14;
		n= 10;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		m = 10;
		n= 4;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		
		m = 273;
		n= 24;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
		
		
		m = 120;
		n= 23;
		print(m,n,extendedGCD(m,n));
		print(m,n,extendedGCD2(m,n));
			
		
	}
}