(数论) 扩展gcd

前言

gcd 可以说是数论中最基础的内容之一。即使不特地学习数论，再学习C语言的时候，或者入门递归的时候必然会学到gcd的求法。

而扩展gcd是在求gcd的同时伴随的求解贝祖定理的两个值，使用场景非常广。如贝祖定理的各种变化，求乘法逆元等操作，是必备的技能。

gcd

扩展gcd前置知识
$$gcd(a,b)=gcd(b,a\%b)$$

一行流模板

int gcd(int a, int b) {
	return b == 0 : a ? gcd(b, a % b);
}
1
2
3

内置函数

// C++ api
#incude <stl_algo.h>

template<typename T>
T __gcd(T __m, T __n);

// C++17
template<typename T>
T gcd(T __m, T __n);
1
2
3
4
5
6
7
8
9

例题

gcd与lcm

$$\begin{gathered} gcd(a,b)=x \\ lcm(a,b)=y \\ a\ast b=gcd(a,b)\ast lcm(a,b)=x\ast y \end{gathered}$$

洛谷：P1029 NOIP2001 普及组 最大公约数和最小公倍数问题

题目描述:

输入两个正整数 $x_0,y_0$，求出满足下列条件的 $P,Q$ 的个数：

1. $P,Q$ 是正整数。

2. 要求 $P,Q$ 以 $x_0$ 为最大公约数，以 $y_0$ 为最小公倍数。

试求：满足条件的所有可能的 $P,Q$ 的个数。

/**
 * https://www.luogu.com.cn/problem/P1029
 * P1029 [NOIP2001 普及组] 最大公约数和最小公倍数问题
 */
#include 
using namespace std;
#define int long long

int gcd(int a, int b) {
    return b == 0 ? a : gcd(b, a % b);
}

signed main() {
    int Gcd, Lcm;
    cin >> Gcd >> Lcm;
    int mult = Gcd * Lcm;

    int ans = 0;
    // 对称性 [1, √mult]
    for (int i = 1; i * i <= mult; i += 1) {
        if (mult % i == 0 && gcd(i, mult / i) == Gcd) {
            ans += 2;
        }
    }
    // 特判，gcd和lcm相等是单个点
    if (Gcd == Lcm) {
        ans += -1;
    }

    cout << ans << endl;
    return 0;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32

扩展gcd

模板

// 注意，这里的x,y是引用
int exgcd(int a, int b, int& x, int& y) {
    if (b == 0) {
        x = 1, y = 0;
        return a;
    }
    int xx, yy;
    int d = exgcd(b, a % b, xx, yy);
    x = yy;
    y = xx - a / b * yy;
    return d;
}
1
2
3
4
5
6
7
8
9
10
11
12

引入与证明

ax + by = gcd(a, b)

核心公式 裴蜀定理_百度百科 (baidu.com)

必存在整数解x,y使得下式成立
$$a\ast x+b\ast y=gcd(a,b)$$

证明

当 $b=0$ 时，$a\ast x+b\ast y=a$ 因此 $x=1y=0(任意)$

当 $b!=0$ 时 根据贝祖定理
{ g c d ( a , b ) = a x + b y ① { g c d ( b , a % b ) = b x 1 + ( a % b ) y 1 ② = b x 1 + ( a − ⌊ a b ⌋ ∗ b ) y 1 = a y 1 + b ( x 1 − a b y 1 ) ① & ② = > { x = y 1 y = x 1 − a b y 1

{gcd(a,b)=ax+by①" role="presentation">{gcd(a,b)=ax+by①$$\{\begin{array}{l}gcd(a,b)=ax+by①\end{array}$$

\\

{gcd(b,a%b)=bx1+(a%b)y1②=bx1+(a−⌊ab⌋∗b)y1=ay1+b(x1−aby1)" role="presentation">⎧⎩⎨⎪⎪gcd(b,a%b)=bx1+(a%b)y1②=bx1+(a−⌊ab⌋∗b)y1=ay1+b(x1−aby1)$$\{\begin{array}{l}gcd(b,a\mathrm{\%}b)=b{x}_1+(a\mathrm{\%}b)y_1②\\ =b{x}_1+(a-\lfloor \frac{a}{b}\rfloor \ast b)y_1\\ =a{y}_1+b(x_1-\frac{a}{b}{y}_1)\end{array}$$

\\ ① \& ②=> \\

{x=y1y=x1−aby1" role="presentation">{x=y1y=x1−aby1$$\{\begin{array}{l}x=y_1\\ y=x_1-\frac{a}{b}{y}_1\end{array}$$

{gcd(a,b)=ax+by①​⎩ ⎨ ⎧​gcd(b,a%b)=bx1​+(a%b)y1​②=bx1​+(a−⌊ba​⌋∗b)y1​=ay1​+b(x1​−ba​y1​)​①&②=>{x=y1​y=x1​−ba​y1​​

因此可以借助基础的gcd运算的同时，把x,y算出

变形

$ a*x + b * y = 1 $

当 等式右边为1 => $a\ast x+b\ast y=1$

即可得a，b两数互质

即求乘法逆元

$ a*x + b * y = c $

当 等式右边为一个常数 => $a\ast x+b\ast y=c$

则相当于在原始两边左右同时乘上了$\frac{c}{gcd(a,b)}$

=> 算出的 $x,y$ 也乘上 $\frac{c}{gcd(a,b)}$ 即可

=> 即求特解

$ a*x + b * y = 0 $

当 等式右边为0 => $a\ast x+b\ast y=0$

首先0也是一个常数，但是一个特殊的常数，因此我们只要凑抵消即可
有特解 { x 0 , y 0 } a ∗ x + b ∗ y = 0 { x = x 0 + b g c d ( a , b ) k y = y 0 − a g c d ( a , b ) k 有特解{\{ x_0, y_0 \} } \\ a*x + b * y = 0 \\

{x=x0+bgcd(a,b)ky=y0−agcd(a,b)k" role="presentation">{x=x0+bgcd(a,b)ky=y0−agcd(a,b)k$$\{\begin{array}{l}x=x_0+\frac{b}{gcd(a,b)}k\\ y=y_0-\frac{a}{gcd(a,b)}k\end{array}$$

有特解{x0​,y0​}a∗x+b∗y=0{x=x0​+gcd(a,b)b​ky=y0​−gcd(a,b)a​k​

=> 即求通解

例题

贝祖定理 (裴蜀定理)

洛谷：P4549 【模板】裴蜀定理

题目描述

给定一个包含 $n$ 个元素的整数序列 $A$，记作 $A_1,A_2,A_3,...,A_n$。

求另一个包含 $n$ 个元素的待定整数序列 $X$，记 $S=\sum _{i=1}^n{A}_i\times X_i$，使得 $S>0$ 且 $S$ 尽可能的小。

---

$ax+by=gcd(a,b)$

扩展，多个整数的叠加 =>

${\sum }_{i=1}^n{a}_i{x}_i=gcd(a_1...a_n)$

/**
 * https://www.luogu.com.cn/problem/P4549
 * P4549 【模板】裴蜀定理
 */
#include 
using namespace std;

int main() {
    int n;
    cin >> n;

    int sum = 0;
    for (int i = 1, a; i <= n; i += 1) {
        cin >> a;
        sum = __gcd(sum, abs(a));
    }
    cout << sum << endl;

    return 0;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20

同余方程

洛谷：P1082 NOIP2012 提高组 同余方程

题目描述

求关于$ x$的同余方程 $ a x \equiv 1 \pmod {b}$ 的最小正整数解。

---

线性同余方程，要么无解，要么有gcd个不用的解

/**
 * https://www.luogu.com.cn/problem/P1082
 * P1082 [NOIP2012 提高组] 同余方程
 * =====================================
 * ax = 1 (mod b)
 * ax + by = 1 求x
 */
#include 
using namespace std;
#define int long long

// 扩展gcd模板
int exgcd(int a, int b, int& x, int& y) {
    if (b == 0) {
        x = 1, y = 0;
        return a;
    }
    int xx, yy;
    int d = exgcd(b, a % b, xx, yy);
    x = yy;
    y = xx - a / b * yy;
    return d;
}

signed main() {
    int a, b;
    cin >> a >> b;

    int x, y;
    exgcd(a, b, x, y);
    // 正整数解
    cout << (x % b + b) % b << endl;

    return 0;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35

费马小定理

$a^{p-1}\equiv 1(\mathrm{mod}p)$ 若p为一个素数，则有

$$\begin{gathered} a^{p-1}\equiv 1(\mathrm{mod}p) \\ a\ast a^{p-2}\equiv 1(\mathrm{mod}p) \end{gathered}$$

这里的x就是a对于mod b的逆元，若b为一个素数，则可以使用费马小定理

$$\begin{gathered} (a/b)\%p= \\ (a\ast b^{-1})\%p= \\ ((a\%p)\ast (b^{p-2}\%p))\%p \end{gathered}$$

练习题：杭电：A/B - 1576

// 费马小定理算除法取模的模板
/**
* molecular 分子
* denominator 分母
*/
int subMod(int molecular, int denominator, int mod) {
    auto binPow = [](int base, int expo, int mod) -> int {
        int ans = 1;
        base %= mod;
        if (expo == 0) {
            ans = 1%mod;
        } else {
            while(expo) {
                if (expo & 1) {
                    ans = ans * base % mod;
                }
                base = base * base % mod;
                expo >>= 1;
            }
        }

        return ans;
    };

    int inverseElement = binPow(denominator, mod-2, mod);
    return (molecular * inverseElement) % mod;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27

不定方程

洛谷：P5656 【模板】二元一次不定方程 (exgcd)

题目描述

给定不定方程

$$ax+by=c$$

若该方程无整数解，输出 $-1$。
若该方程有整数解，且有正整数解，则输出其正整数解的数量，所有正整数解中 $x$ 的最小值，所有正整数解中 $y$ 的最小值，所有正整数解中 $x$ 的最大值，以及所有正整数解中 $y$ 的最大值。
若方程有整数解，但没有正整数解，你需要输出所有整数解中 $x$ 的最小正整数值， $y$ 的最小正整数值。

正整数解即为 $x,y$ 均为正整数的解，$0$ 不是正整数。
整数解即为 $x,y$ 均为整数的解。
$x$ 的最小正整数值即所有 $x$ 为正整数的整数解中 $x$ 的最小值，$y$ 同理。

---

本题非常绕，流程就是先算特解，再算通解，再将x化为最小正整数，再判断y

推荐博客：洛谷P5656 【模板】二元一次不定方程 (exgcd) 题解_q779的博客-CSDN博客

/**
 * https://www.luogu.com.cn/problem/P5656
 * P5656 【模板】二元一次不定方程 (exgcd)
 */
#include 
using namespace std;
#define int long long

int exgcd(int a, int b, int &x, int &y) {
    if (b == 0) {
        x = 1, y = 0;
        return a;
    }
    int xx, yy;
    int d = exgcd(b, a % b, xx, yy);
    x = yy;
    y = xx - (a / b) * yy;
    return d;
}

signed main() {
    int n;
    scanf("%lld", &n);

    for (int i = 1, a, b, c, x, y; i <= n; i += 1) {
        scanf("%lld %lld %lld", &a, &b, &c);
        int GCD = exgcd(a, b, x, y);

        // 无解
        if (c % GCD) {
            printf("-1\n");
            continue;
        }
        // exgcd 求的是 ax + by = gcd
        // 现在是 ax + by = c

        // 化特解
        // 等式右边gcd => c
        // x,y化为特解
        x *= (c / GCD);
        y *= (c / GCD);

        // 化通解
        // x = x0 + b/gcd * k
        // y = y0 - a/gcd * k
        int bgcd = b / GCD;
        int agcd = a / GCD;
        // 找最小正整数k，目的是让x为正
        // 使得 x0 + b/gcd * k >= 1
        int k = ceil((1.0 - x) / bgcd);
        // 化为通解
        x = x + bgcd * k;
        y = y - agcd * k;

        // 此时x已为正，仅需考虑y
        // 若 y <= 0 则仅存在整数解
        if (y <= 0) {
            // 用上面求x的方式再求一下y
            k = ceil((1.0 - y) / agcd);
            y = y + agcd * k;
            // x,y 的最小正整数解
            printf("%lld %lld\n", x, y);
        }
        // x,y 均为正整数
        // 此时x为最小正整数，因此用y计算
        else {
            // 正整数解的个数
            // 向上取整
            printf("%lld ", (y - 1) / agcd + 1);
            // x 的最小正整数值
            // 前面确保了x最小
            printf("%lld ", x);
            // y 的最小正整数值
            // -1%+1
            printf("%lld ", (y - 1) % agcd + 1);
            // x 的最大正整数值
            // 累计组数(y)*组值(x)
            printf("%lld ", x + (y - 1) / agcd * bgcd);
            // y 的最大正整数值
            // 前面确保了y最小
            printf("%lld\n", y);
        }
    }
    return 0;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85

线性求逆元

洛谷：P3811 【模板】乘法逆元

题目描述

给定 $n,p$ 求 $1\sim n$ 中所有整数在模 $p$ 意义下的乘法逆元。

这里 $a$ 模 $p$ 的乘法逆元定义为 $ax\equiv 1(\mathrm{mod}p)$ 的解。

---

线性求逆元并没有用到扩展gcd，但是往往是一个由扩展gcd拓展出来的问题（就像费马小定理一样）

证明：设求第i个数的逆元
设 p ≡ k ∗ i + r ( p 为素数， k 为倍数， i 为要求第几个数的逆元， r 为余数 ) 变形 k ∗ i + r ≡ 0 ( m o d p ) 两边同时乘 i − 1 和 r − 1 构造出 i − 1 k ∗ r − 1 + i − 1 ≡ 0 ( m o d p ) 移项 i − 1 ≡ − k ∗ r − 1 ( m o d p ) 用 p 和 i 表示倍数 k 和余数 r { k = p / i r = p % i i − 1 ≡ − ⌊ p / i ⌋ ∗ ( p % i ) − 1 ( m o d p ) 设 \quad p \equiv k * i + r \quad (p为素数，k为倍数，i为要求第几个数的逆元，r为余数) \\ 变形 \quad k * i + r \equiv 0 \pmod p \\ 两边同时乘 i^{-1} 和 r^{-1} 构造出i^{-1} \\ k * r^{-1} + i^{-1} \equiv 0 \pmod p \\ 移项 \quad i^{-1} \equiv - k * r^{-1} \pmod p \\ 用p和i表示倍数k和余数r

{k=p/ir=p%i" role="presentation">{k=p/ir=p%i$$\{\begin{array}{l}k=p/i\\ r=p\mathrm{\%}i\end{array}$$

\\ i^{-1} \equiv - \lfloor p / i \rfloor * (p \% i)^{-1} \pmod p 设p≡k∗i+r(p为素数，k为倍数，i为要求第几个数的逆元，r为余数)变形k∗i+r≡0(modp)两边同时乘i−1和r−1构造出i−1k∗r−1+i−1≡0(modp)移项i−1≡−k∗r−1(modp)用p和i表示倍数k和余数r{k=p/ir=p%i​i−1≡−⌊p/i⌋∗(p%i)−1(modp)

/**
 * https://www.luogu.com.cn/problem/P3811
 * P3811 【模板】乘法逆元
 * 线性求逆元
 */
#include 
using namespace std;
#define int long long

signed main() {
    int n, p;
    scanf("%lld %lld", &n, &p);

    // i^-1 = -floor(p/i) * (p%i)^-1 % p
    vector<int> inv(n + 1);
    inv[1] = 1;
    for (int i = 2; i <= n; i += 1) {
        // p/i化为正数
        inv[i] = (p - (p / i)) * inv[p % i] % p;
    }

    for (int i = 1; i <= n; i += 1) {
        printf("%lld\n", inv[i]);
    }

    return 0;
}
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27

---

---

---

END