GCD+拓展欧几里德

GCD

求Gcd

两种常用方法，均基于辗转相除，第二种方法在某些问题的处理上可能出锅。

inline int gcd(int x,int y){
    return y==0?x:gcd(y,x%y);
}
inline int _gcd(int x,int y){
    while(y^=x^=y^=x%=y);
	return x;
}

求LCM

$gcd(a,b)lcm(a,b)=ab$

拓展欧几里德(exgcd)

Native

接触拓展欧几里德之前，数论接触的并不多，尽量把式子展开看好了。

首先有两个定理：

裴蜀定理：对于方程 $ax+by=c \ \ \ x,y\in \mathbb{Z}^+$，当且仅当 $gcd(a,b)|c$ 时，有整数解。

欧几里德：$gcd(a,b)=gcd(b,a \ mod \ b)$。

所以首先对于方程 $ax+by=gcd(a,b)$，可以化为：

$$ax+by=gcd(b,a \ mod \ b)$$ $$bx^*+(a-\lfloor\frac{a}{b}\rfloor b)y^*=gcd(b,a \ mod \ b)$$

因此带换回去：

$$ax+by=bx^*+(x^*-\lfloor\frac{a}{b}\rfloor y^*)y^*$$

一一对应得到

$$x=y^*,y=(x^*-\lfloor\frac{a}{b}\rfloor y^*)$$

因此可以递归求解。

再看方程 $ax+by=c$，如果有解，那么：

$$\frac{ax}{gcd(a,b)}+\frac{by}{gcd(a,b)}=\frac{gcd(a,b)}{gcd(a,b)}$$ $$ax\frac{c}{gcd(a,b)}+by\frac{c}{gcd(a,b)}=c$$

因此求出此时的 $x,y$ 再乘以 $\frac{c}{gcd(a,b)}$ 就可以得到该方程的解。

解的调整

注意到拓展欧几里得算法只求出了$x$的最小整数解，没有求出完整的解集。

我们可以对解进行调整：

$$ax+abk+by-abk=c$$ $$a(x+bk)+b(y-ak)=c$$

同理可得：

$$a(x-bk)+b(y+ak)=c$$

因此调整公式为：

$$a(x\pm bk)+b(y\mp ak)=c,k\in\mathbb{Q}$$

问题在于如何选取 $k$ 保证 $x,y$ 为整数，我们需要保证的是 $bk,ak\in \mathbb{Z}$，令 $bk,ak$ 取得最小的正整数（一般讨论的方程中 $a,b$ 均为正，若不为正此处需分类讨论），显然发现 $k=\frac{1}{(a,b)}$ 时有 $bk_{min}=\frac{b}{(a,b)},ak_{min}=\frac{a}{(a,b)}$（可作为定理记忆，此处对证明不做多赘述），而我们调整解时需满足 $bk_{min}|bk$ 与 $ak_{min}|ak$，以此即可得出通解形式：

$$x*=x\pm s\frac{b}{(a,b)},y*=y\mp s\frac{a}{(a,b)},s\in\mathbb{Z}$$

接下来以此通解形式调整解即可得出 $x,y$ 对应的最大整数解和最小整数解。

同时可以注意到方程的正整数解个数计算方式是 $\frac{x_{max}-x_{min}}{b}+1$，其中 $x_{max}$ 为最大整数解，而 $x_{min}$ 为最小整数解，该公式可以理解为一个已知两端公差为 $b$ 的等差数列的项数计算。

Code

拓展欧几里德：

void exgcd(int a,int b,int &x,int &y){
	if(!b){x=1,y=0;return;}
	exgcd(b,a%b,x,y);
	int temp=x;x=y;y=temp-(a/b)*y;
}
int exgcd(int a,int b,int &x,int &y){//求gcd
	if(!b){x=1,y=0;return a;}
	int gcd=exgcd(b,a%b,x,y);
	int temp=x;x=y;y=temp-(a/b)*y;
	return gcd;
}

拓展欧几里德解线性同余方程：

void exgcd(int a,int b,int &x,int &y){
	if(!b){x=1,y=0;return;}
	exgcd(b,a%b,x,y);
	int temp=x;x=y;y=temp-(a/b)*y;
}
bool solve(int a,int b,int c,int &x,int &y){
	int d=gcd(a,b);
	if(c%d)return false;
	exgcd(a,b,x,y);
	int k=c/d;
	x*=k,y*=k;
	return true;
}

或：

int exgcd(int a,int b,int &x,int &y){
	if(!b){x=1,y=0;return a;}
	int gcd=exgcd(b,a%b,x,y);
	int temp=x;x=y;y=temp-(a/b)*y;
	return gcd;
}
bool solve(int a,int b,int c,int &x,int &y){
	int d=exgcd(a,b,x,y);
	if(c%d)return false;
	int k=c/d;
	x*=k,y*=k;
	return true;
}