经典算法（5）- 用二进制方法实现扩展的最大公约数(Extended GCD)

二进制方法中，只需要移位（<<和>>）和加减操作（+和-），不像欧几里德算法中需要乘法和除法运算。虽然算法效率更高，但是程序的可读性和可维护性差一些。

如果设d=gcd(u,v) = u.x + v.y, 本算法涉及到六种操作：

1）已知ext_gcd(u,v)如何求ext_gcd(u,2v)=u'.x' + v'.y'，其中u为奇数，v可奇可偶，d=gcd(u,v)为奇数；

2）已知ext_gcd(u,v)如何求ext_gcd(2u,v)=u'.x' + v'.y'，其中v为奇数，u可奇可偶，d=gcd(u,v)为奇数；

3）已知ext_gcd(u-v,v)如何求ext_gcd(u,v)=u'.x' + v'.y'；

4）已知ext_gcd(u,u-v)如何求ext_gcd(u,v)=u'.x' + v'.y'；

5）已知u/2^c = v = d（c为大于0的整数），如何求ext_gcd(u,v)=u.x' + v.y'；

6）已知u= v/2^c = d（c为大于0的整数），如何求ext_gcd(u,v)=u.x' + v.y'.

其中第1）种操作比较麻烦，第2）种操作只需要在第1）种操作的基础上将u和v交换一下。下面介绍一下第1）种操作的原理：

在第1）操作中，因为d=gcd(u,2v)=gcd(u-v,v)=(u-v)x+v(x+y)，设u1=u-v，x1=x，v1=v，y1=x+y，即问题转化为已知d=u1.x1+u1.y1，求ext_gcd(u,2v)=u.x' + (2v)y'=d；

如果y为偶数，显然可以采用这样一组值：x'=x，y'=y/2；但是如果y为奇数，则需要将d表示为d=u1(x1+v1.k) + v1(y1-u1.k)，其中k为正奇数1，3，5....，找到这样的k，满足x' = x +v.k和y'=(y-u.k)/2都不为0。

第3）种操作比较简单，因为d=gcd(u-v,v) = (u-v)x + vt = ux + v(y-x)，即x'=x，y'=y-x;

第4）种操作类似，因为d=gcd(u,u-v) = ux + (u-v)t = u(x+y) - vy，即x'=x+y，y'=-y;

第5）种操作,可以得到d=v=u + (1-2^c).v，即x'=1，y'=1-2^c；

第6）种操作只需要在第5）种操作的基础上将u和v交换一下。

在下面的实现中，只有上面第2）和6）种操作需要改变v和y的符号。

import java.util.Stack;
/**
*
* @author ljs 2011-5-19
*
* solve extended gcd using Binary Method
*
*/
public class ExtGCD_Binary {
/**
* caculate x,y to satisfy Bezout Identity: u.x + v.y = gcd(u,v)
* @param u
* @param v
* @return {d,x,y} for: d = u.x + v.y;
*/
public static int[] extGCDBinary(int u,int v){
u=(u<0)?-u:u;
v=(v<0)?-v:v;

if(u==0)
return new int[]{v,0,1};
if(v==0)
return new int[]{u,1,0};

int k=0;
while((u & 0x01)==0 && (v & 0x01) == 0){
u>>=1; //divide by 2
v>>=1;
k++;
}
//with respect to d=gcd(u,v)=u.s+v.t
/* the operation type is defined as:
1: c, u / 2^c (push two numbers to stack: first push c, then push the op number 1)
2: c, v / 2^c  (push two numbers to stack: first push c, then push the op number 2)
3: gcd(u-v, v) (only push the op number 3 to stack)
4: gcd(u, u-v) (only push the op number 4 to stack)
*/
Stack<Integer> ops= new Stack<Integer>();
//at this time, there is at least one number is odd between m and n
int t=-v; //set it negative for later comparison of (t>0)
if((v & 0x01)==1){
//if v is odd
t = u;
}

//process t as a possible even number
boolean firstNoOp = true;
while(t != 0){
if(firstNoOp) {
firstNoOp = false;
}else{
if(t>0)
ops.push(3);
else
ops.push(4);
}

int c=0;
while((t & 0x01)==0){
//do until t is not even
t>>=1;
c++;
}
if(t>0) {//u > v (the max is replaced by |t|)
u=t;
if(c>0) {
ops.push(c);
ops.push(1);
}
}else{ //u<v (the max is replaced by |t|)
v=-t;
if(c>0){
ops.push(c);
ops.push(2);
}
}
//now u and v are all odd, then u-v is even
t = u-v;
}
int d = u;

//the following steps are to caculate x,y in d=u.x + v.y

int x=1,y=0;
//special processing for the last operation
if(!ops.isEmpty()){
//e.g. for input u=3,v=3, then ops is empty, d=3
int op = ops.pop();
int c = ops.pop();
//op number 3 and 4 can not be the last step since: u-v = v, then v is even; or u=u-v, then u=0;
//both are not possible
int tmp = 1<<c;
if(op == 1){
x = 1;
y = 1- tmp;

u = v;
for(int i=0;i<c;i++)
u <<= 1;
}else if(op == 2){
x = 1- tmp;
y = -1;

v = -u;
for(int i=0;i<c;i++)
v <<= 1;
}
}
while(!ops.isEmpty()){
int op = ops.pop();
switch(op){
case 1:{
int c = ops.pop();
for(int i=0;i<c;i++){
if((x & 0x01)==0){
//if x is even
x >>= 1;
}else{
y += u;
int tmp = x - v;
while(y ==0 || tmp ==0){
y += (u<<1); //*2
tmp -= (v<<1); //*2
};
x = tmp >>1;
}
u <<= 1;
}
break;
}
case 2:{
int c = ops.pop();
for(int i=0;i<c;i++){
if((y & 0x01)==0){
//if y is even
y >>= 1;
}else{
x += v;
int tmp = y-u;
while(x ==0 || tmp ==0){
x += (v<<1);  //*2
tmp -= (u<<1);  //*2
};
y = tmp >>1;
}
v <<= 1;
}
y = -y;
v = -v;
break;
}
case 3:
y -= x;
u += v;
break;
case 4:
x += y;
v = u - v;

y = -y;
break;
}
}
//e.g. gcd(1,4)
y = (v<0)?-y:y;

for(int i=0;i<k;i++)
d <<=1;
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,extGCDBinary(m,n));

//co-prime
m = 15;
n= 28;
print(m,n,extGCDBinary(m,n));

m = 6;
n= 3;
print(m,n,extGCDBinary(m,n));

m = 6;
n= 3;
print(m,n,extGCDBinary(m,n));

m = 6;
n= 0;
print(m,n,extGCDBinary(m,n));

m = 0;
n= 6;
print(m,n,extGCDBinary(m,n));

m = 0;
n= 0;
print(m,n,extGCDBinary(m,n));

m = 1;
n= 1;
print(m,n,extGCDBinary(m,n));

m = 3;
n= 3;
print(m,n,extGCDBinary(m,n));

m = 2;
n= 2;
print(m,n,extGCDBinary(m,n));

m = 1;
n= 4;
print(m,n,extGCDBinary(m,n));

m = 4;
n= 1;
print(m,n,extGCDBinary(m,n));

m = 10;
n= 14;
print(m,n,extGCDBinary(m,n));

m = 14;
n= 10;
print(m,n,extGCDBinary(m,n));

m = 10;
n= 4;
print(m,n,extGCDBinary(m,n));

m = 273;
n= 24;
print(m,n,extGCDBinary(m,n));

m = 120;
n= 23;
print(m,n,extGCDBinary(m,n));

}
}