欧几里得除法以及拓展

最大公约数

暴力

时间复杂度$O(\min (a, b))$

辗转相除

辗转相除法，又名欧几里得除法，目的是求出两个正整数的最大公约数。

两个正整数a和b（a>b），它们的最大公约数等于a除以b的余数c和b之间的最大公约数。比如10和25，25除以10商2余5,那么10和25的最大公约数，等同于10和5的最大公约数。

int gcd(int a, int b) {
  if (!b || a == b) return a;
  if (a > b) return gcd(b, a % b);
  else return gcd(a, b % a);
}

时间复杂度近似$O(\log \min(a, b))$

更相减损术

辗转相除法涉及取模运算，当整数较大时，取模运算的效率比较差。这时可以看更相减损术。

更相减损术， 出自于中国古代的《九章算术》，也是一种求最大公约数的算法。两个正整数a和b（a>b），它们的最大公约数等于a-b的差值c和较小数b的最大公约数。比如10和25，25减去10的差是15,那么10和25的最大公约数，等同于10和15的最大公约数。

int gcd(int a, int b) {
  if (a == b) return a;
  if (a > b) return gcd(b, a - b);
  else return gcd(a, b - a);
}

时间复杂度$O(\max (a, b)))$

辗转相除结合更相减损术

更相减损术不稳定，运算次数更大，尤其是当两个数比较悬殊的时候。考虑两者结合：

a偶b偶，$gcd(a,b) = 2*gcd(a/2, b/2) = 2*gcd(a>>1, b>>1)$

a偶b奇，$gcd(a,b) = gcd(a/2, b) = gcd(a>>1, b)$

a奇b偶，$gcd(a,b) = gcd(a, b/2) = gcd(a, b>>1)$

a奇b奇，利用更相减损术运算一次，$gcd(a,b) = gcd(b, a-b)$， 此时a-b必然是偶数，又可以继续进行移位运算。

int gcd(int a, int b) {
  if (a == b) return a;
  if (a & 1) {
    if (b && 1) return gcd(a, std::abs(a - b));
    else return gcd(a, b >> 1);
  }
 	else {
    if (b & 1) return gcd(a >> 1, b);
    else return 2 * gcd(a >> 1, b >> 1);
  }
}

时间复杂度$O(\log \max (a, b)))$

参考

https://houbb.github.io/2017/08/23/math-03-common-gcd-03

https://zhuanlan.zhihu.com/p/31824895

最小公倍数

a * b / gcd(a, b)

扩展欧几里得

原理

可用于求 $a x + b y = gcd(a,\ b)$ 的一组解$(x_1,\ y_1)$。

$$a x_1 + b y_1 = gcd(a,\ b) \tag{1}$$

下面是简要求解过程。

用b与$a\ \%\ b$来替代a和b，得到

$$bx_2 + (a\ \%\ b)\ y_2 = gcd (b,\ a\ \%\ b) \tag{2}$$

根据欧几里得算法以及$a\ \%\ b$的算数表达形式

$$\begin{align} gcd(a,\ b)\ &=\ gcd (b,\ a\ \%\ b) \tag{3} \\ a\ \%\ b\ &=\ a\ -\ \lfloor \frac{a}{b} \rfloor b \tag{4} \end{align}$$

联立式子$(2)$、$(3)$和$(4)$

$$a y_2\ +\ b \left [ x_2\ -\ \lfloor \frac{a}{b} \rfloor y_2 \right ]\ =\ gcd (a, b) \tag{5}$$

联立式子$(1)$和$(5)$

$$\begin{align} x_1\ &=\ y_2 \\ y_1\ &=\ x_2\ -\ \lfloor \frac{a}{b} \rfloor y_2 \end{align}$$

由此写出以下程序。

int exgcd(int a,int b,int &x,int &y)
{
    if (b==0) {
        x = 1;
        y = 0;
        return a;
    }
    int r = exgcd(b, a%b, x, y);
    int t = y;
    y = x - (a / b) * y;     //2的情况
    x = t;
    return r;
}

通解

如果 $(x_0, y_0)$ 是 $a x + b y = c$ 一组解，那么通解可以表示为

$$\begin{align} x\ &=\ x_0\ +\ \frac{b}{c} k \\ y\ &=\ y_0\ -\ \frac{a}{c} k \end{align}$$

事实上

$$a x + b y\ =\ a \left( x_0\ +\ \frac{b}{c} k \right)\ +\ b \left( y_0\ -\ \frac{a}{c} k \right)\ =\ a x_0 + b y_0\ =\ c$$

参考

https://zhuanlan.zhihu.com/p/100567253

https://zhuanlan.zhihu.com/p/42707457

https://www.desgard.com/algo/docs/part2/ch02/3-ext-euclidean/

https://ksmeow.moe/euclid/

https://www.cnblogs.com/wkfvawl/p/9350867.html

例子

以Acwing 203 同余方程为例，$ax \equiv 1\ (\mod m)$ 可以写成 $a x - 1 = m y$，即 $ax - my = 1$。形式上等同于 $a x + b y = c,\ b = -m,\ c = 1$。

根据扩展欧几里得算法可以解出一个$x_0$，那么根据通项解可知，就是在$x_0$上加减若干个 $\frac{b}{c}m = m$ 。由于$x_0$可能是负数，或者是大于$m$的正数，通过$(x \% b + b) \% b$使其成为$(0, m)$之间的正整数。

#include <cstdio>
using namespace std;

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), temp = y;
    y = x - a / b * y;
    x = temp;
    return gcd;
}
int main() {
    int a, b, x, y;
    scanf("%d%d", &a, &b);
    exgcd(a, b, x, y);
    printf("%d\n", (x % b + b) % b);
    return 0;
}

其他还有

• Leetcode365 水壶问题