算法实现/数学/扩展欧几里得算法
外观
#include <stdio.h>
#include <stdlib.h> // strtol() function
void xgcd(long *result, long a, long b){
long aa[2]={1,0}, bb[2]={0,1}, q;
while(1) {
q = a / b; a = a % b;
aa[0] = aa[0] - q*aa[1]; bb[0] = bb[0] - q*bb[1];
if (a == 0) {
result[0] = b; result[1] = aa[1]; result[2] = bb[1];
return;
};
q = b / a; b = b % a;
aa[1] = aa[1] - q*aa[0]; bb[1] = bb[1] - q*bb[0];
if (b == 0) {
result[0] = a; result[1] = aa[0]; result[2] = bb[0];
return;
};
};
}
int main(int argc, char** argv){
long a, b, c[3];
char *end;
if(argc < 3) {printf("Usage: %s a b (compute xgcd(a,b))\n",argv[0]); exit(-1);};
printf("long int is %lu bytes (or %lu bits)\n",sizeof(long),8*sizeof(long));
a = strtol(argv[1], &end, 10);
b = strtol(argv[2], &end, 10);
printf("xgcd(%ld,%ld) = ",a,b);
xgcd(c,a,b);
printf("%ld = %ld*%ld + %ld*%ld\n",c[0],c[1],a,c[2],b);
return 0;
}
//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;
}
public class xgcd {
public static long[] xgcd(long a, long b){
long[] retvals={0,0,0};
long aa[]={1,0}, bb[]={0,1}, q=0;
while(true) {
q = a / b; a = a % b;
aa[0] = aa[0] - q*aa[1]; bb[0] = bb[0] - q*bb[1];
if (a == 0) {
retvals[0] = b; retvals[1] = aa[1]; retvals[2] = bb[1];
return retvals;
};
q = b / a; b = b % a;
aa[1] = aa[1] - q*aa[0]; bb[1] = bb[1] - q*bb[0];
if (b == 0) {
retvals[0] = a; retvals[1] = aa[0]; retvals[2] = bb[0];
return retvals;
};
}
}
public static void main(String[] args){
if(args.length < 2){
System.out.println("Usage: xgcd(a,b)"); System.exit(-1);};
long a = Integer.parseInt(args[0]);
long b = Integer.parseInt(args[1]);
long[] retvals;
retvals = xgcd(a,b);
System.out.println("xgcd("+args[0]+", "+args[1]+") = ["
+String.valueOf(retvals[0])+", "+
String.valueOf(retvals[1])+", "+
String.valueOf(retvals[2])+"]");
System.out.println(" "+String.valueOf(retvals[0])+" = "
+args[0]+"*("+String.valueOf(retvals[1])+") + "
+args[1]+"*("+String.valueOf(retvals[2])+")");
};
}
import java.math.BigInteger;
public class XGCD {
// Assumes a and b are positive
public static BigInteger[] xgcd(BigInteger a, BigInteger b) {
BigInteger x = a, y=b;
BigInteger[] qrem = new BigInteger[2];
BigInteger[] result = new BigInteger[3];
BigInteger x0 = BigInteger.ONE, x1 = BigInteger.ZERO;
BigInteger y0 = BigInteger.ZERO, y1 = BigInteger.ONE;
while (true) {
qrem = x.divideAndRemainder(y);
x = qrem[1];
x0 = x0.subtract(y0.multiply(qrem[0]));
x1 = x1.subtract(y1.multiply(qrem[0]));
if (x.equals(BigInteger.ZERO)) {
result[0] = y;
result[1] = y0;
result[2] = y1;
return result;
}
qrem = y.divideAndRemainder(x);
y = qrem[1];
y0 = y0.subtract(x0.multiply(qrem[0]));
y1 = y1.subtract(x1.multiply(qrem[0]));
if (y.equals(BigInteger.ZERO)) {
result[0] = x;
result[1] = x0;
result[2] = x1;
return result;
}
}
}
public static void main(String[] args) {
BigInteger a = new BigInteger(args[0]);
BigInteger b = new BigInteger(args[1]);
BigInteger[] result = new BigInteger[3];
System.out.println("a = " + a.toString() + "; b = " + b.toString());
result = xgcd(a,b);
System.out.println("xgcd = " + result[0].toString() +
" [" + result[1].toString() + ", " + result[2].toString() + "]");
}
}
这两个函数都以正整数 a、b 作为输入,并返回三元组 (g, x, y),使得 ax + by = g = gcd(a, b)。
from typing import Tuple
def xgcd(a: int, b: int) -> Tuple[int, int, int]:
"""return (g, x, y) such that a*x + b*y = g = gcd(a, b)"""
x0, x1, y0, y1 = 0, 1, 1, 0
while a != 0:
(q, a), b = divmod(b, a), a
y0, y1 = y1, y0 - q * y1
x0, x1 = x1, x0 - q * x1
return b, x0, y0
import sys
from typing import Tuple
sys.setrecursionlimit(1000000) # long type, 32bit OS 4B, 64bit OS 8B (1bit for sign)
def egcd(a: int, b: int) -> Tuple[int, int, int]:
"""return (g, x, y) such that a*x + b*y = g = gcd(a, b)"""
if a == 0:
return (b, 0, 1)
else:
b_div_a, b_mod_a = divmod(b, a)
g, x, y = egcd(b_mod_a, a)
return (g, y - b_div_a * x, x)
扩展 GCD 算法在寻找 模逆 中的应用
def modinv(a: int, b: int) -> int:
"""return x such that (x * a) % b == 1"""
g, x, _ = egcd(a, b)
if g != 1:
raise Exception('gcd(a, b) != 1')
return x % b
(defn xgcd
"Extended Euclidean Algorithm. Returns [gcd(a,b) x y] where ax + by = gcd(a,b)."
[a b]
(if (= a 0)
[b 0 1]
(let [[g x y] (xgcd (mod b a) a)]
[g (- y (* (Math/floorDiv b a) x)) x])))
euclid :: Integral a => a -> a -> a
euclid 0 b = abs b
euclid a 0 = abs a
euclid a b = euclid b $ rem a b
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)