数论——同余问题全家桶2 不定方程和同余方程

不定方程

不定方程又称丢番图方程。这里只讨论二元一次不定方程$ax+by=c$，也就是算法竞赛常考的那部分。

不定方程$ax+by=c$中$x$和$y$都是未知数，因此这个方程又可以转换成$ax\equiv c(\mathrm{mod}b)$。

裴蜀定理

裴蜀定理又称贝祖定理：

• 一定存在整数 $x,y$，满足 $ax+by=\text{gcd}(a,b)$。

可以用来判定不定方程是否存在整数解。

例如 $6x+8y=\text{gcd}(6,8)=2$，一定有整数解。

其中 $x=-1,y=1$ 就是一组解。

【裴蜀定理的推论】

1. 一定存在整数 $x,y$，满足 $ax+by=\text{gcd}(a,b)\times n$

因此，对于一个二元一次不定方程 $ax+by=c$，当 $\text{gcd}(a,b)|c$ 时，有解。

2. 一定存在整数 $x_1,x_2,x_3,\dots$，满足 $a_1{x}_1+a_2{x}_2+a_3{x}_3+\dots a_k{x}_k=\text{gcd}(a_1,a_2,a_3,\dots a_k)\times n$。

$a,b$ 的正负是不影响结果的。因为 $a,b$ 如果存在解，那么 $a,-b$ 也一定存在解，只不过在 $a,b$ 解的基础上添上一个负号。

因此，在用欧几里得算法求 $\text{gcd}(a,b)$ 时，为了避免出现负数，可以带入 $a,b$ 的绝对值。

裴蜀定理OJ示例

这个定理只给出了不定方程的解的存在，所以这个OJ考的其实是其他知识。

【模板】裴蜀定理 - 洛谷

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

根据推论$a_1{x}_1+a_2{x}_2+a_3{x}_3+\dots a_k{x}_k=\text{gcd}(a_1,a_2,a_3,\dots a_k)\times n$，方程的最小值就是求所有数的$\text{gcd}$。

系数这里选择全改成正数。

#include<bits/stdc++.h>
using namespace std;int gcd(int a, int b) {return b ? gcd(b, a % b) : a;
}void ac() {int n; cin >> n;int ans = 1, x;--n; cin >> ans;while (n--) {cin >> x;ans = gcd(abs(ans), abs(x));}cout << ans << endl;
}int main() {int T = 1;//cin >> T;while (T--)ac();return 0;
}

扩展欧几里得算法

对$ax+by=\text{gcd}(a,b)$，当 $b=0$ 时，$ax+by=a$，因此有一组解为：$x=1,y=0$；

当 $b\ne 0$ 时，

• 由欧几里得算法得：$\text{gcd}(a,b)=\text{gcd}(b,a\mathrm{mod}b)$

• 由裴蜀定理得：
  $\text{gcd}(a,b)=ax+by$

  $\text{gcd}(b,a\mathrm{mod}b)=b{x}_1+(a\mathrm{mod}b)y_1$
  $=b{x}_1+(a-\frac{a}{b}\times b)y_1$（在整数除法中存在舍弃小数部分的行为，不能约掉分母）
  $=a{y}_1+b(x_1-\frac{a}{b}{y}_1)$

因为$\text{gcd}(a,b)=\text{gcd}(b,a\mathrm{mod}b)$，所以：
$$\{\begin{array}{l}x=y_1,\\ y=x_1-\frac{a}{b}{y}_1\end{array}$$

这组解需要根据递归和回溯过程来理解用法。

扩展欧几里得算法示例：

//可改成其他类型，用模板是为了方便拷贝
template<class intT>
intT exgcd(intT a, intT b, intT& x, intT& y) {if (b == 0) {x = 1, y = 0;return a;}intT x1, y1, d;d = exgcd(b, a % b, x1, y1);x = y1, y = x1 - a / b * y1;return d;
}

例如：$\text{gcd}(6,8,x,y)\to \text{gcd}(8,6,x,y)\to \text{gcd}(6,2,x,y)\to \text{gcd}(2,0,1,0)$，

回溯：$\text{gcd}(2,0,1,0)\to \text{gcd}(6,2,0,1)\to \text{gcd}(8,6,1,-1)\to \text{gcd}(6,8,-1,1)$。

所以不定方程$6x+8y=\text{gcd}(6,8)$的1个特解是$x=-1,y=1$。

于是可以利用递归，先求出下一层的 $x_1,y_1$，再求出当前层的 $x,y$。

上述递归过程，可以求出一组特解：$x_0,y_0$。

然后构造出通解为：
$\{\begin{array}{l}x=x_0+\frac{b}{\text{gcd}(a,b)}\times k\\ y=y_0-\frac{a}{\text{gcd}(a,b)}\times k\end{array}$
其中 $k$ 是正整数。

通解的构造：$a(x_0+\Delta x)+b(y_0-\Delta y)=\text{gcd}(a,b)$，

将左边拆开：$a{x}_0+b{y}_0+a\Delta x-b\Delta y=\text{gcd}(a,b)$，

为了和原等式对应，$a\Delta x-b\Delta y=0$，所以$\frac{\Delta x}{\Delta y}=\frac{b}{a}$，

将$\frac{b}{a}$化简得$\frac{\Delta x}{\Delta y}=\frac{\frac{b}{\text{gcd}(a,b)}}{\frac{a}{\text{gcd}(a,b)}}$。所以$\Delta x$、$\Delta y$可看成$\frac{b}{\text{gcd}(a,b)}$和$\frac{a}{\text{gcd}(a,b)}$。

为了得到无穷的通解，引入正整数$k$，

所以就有通解：
$\{\begin{array}{l}x=x_0+\frac{b}{\text{gcd}(a,b)}\times k\\ y=y_0-\frac{a}{\text{gcd}(a,b)}\times k\end{array}$

时间复杂度：与欧几里得算法得时间复杂度一致，为 $O(\log n)$。

扩欧求二元一次不定方程的解

根据扩展欧几里得算法，求解二元一次不定方程 $ax+by=c$ 的解的流程：

1. 利用扩欧算法求出不定方程 $ax+by=\gcd (a,b)$ 一组特解 $x_0,y_0$，以及系数的最大公约数 $d$；

2. 由裴蜀定理判断是否有解：
   a. 如果 $d\nmid c$（$d$不是$c$的约数）：无解；
   b. 如果 $d\mid c$：有解。

3. 如果有解：
   a. 原方程的一组特解为：$x_1=\frac{x_0\times c}{d},y_1=\frac{y_0\times c}{d}$；

   根据裴蜀定理的推论，$c$可能不等于$\text{gcd}(a,b)$，所以$a{x}_0+b{y}_0=\text{gcd}(a,b)$，记$d=\text{gcd}(a,b)$，

   则$\frac{c}{d}(a{x}_0+b{y}_0)=d\times \frac{c}{d}$，

   所以$a(\frac{x_0\times c}{d})+b(\frac{y_0\times c}{d})=d\times \frac{c}{d}$，

   因为$d|c$，所以$a(\frac{x_0\times c}{d})+b(\frac{y_0\times c}{d})=c$，

   即方程$ax+by=c$的一组特解为：$x_1=\frac{x_0\times c}{d},y_1=\frac{y_0\times c}{d}$。

   b. 根据特解可得通解：$x=x_1+\frac{b}{\gcd (a,b)}\times k,y=y_1-\frac{a}{\gcd (a,b)}\times k$。

扩欧定理OJ示例

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

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

借助这题将不定方程的解尽可能多的解释细节。

解决这个题需要使用扩欧算法解二元一次不定方程：

1. 先用exgcd求出$ax+by=\text{gcd}(a,b)$的一组特解，记$d=\text{gcd}(a,b)$。

2. $d\nmid c$：无解，否则有解。

3. 根据上文推论，方程$ax+by=c$的一组特解为$x_1=\frac{x_0\times c}{d},y_1=\frac{y_0\times c}{d}$。利用特解可以求得通解$x=x_1+\frac{b}{\gcd (a,b)}\times k,y=y_1-\frac{a}{\gcd (a,b)}\times k$。

   • 观察通解可发现，当$k$逐渐增加，$x$也逐渐增加，而$y$反而减小，反之同理。
     所以先将$x$补正变成最小的正整数，此时$y$即为最大值。

     • $y$为最大值时，当$y\le 0$时，方程无正整数解。
     • $y$为最大值时，当$y>0$时，方程有正整数解。

     关于如何补正：$x=x_1+\frac{b}{\gcd (a,b)}\times k$和$5=1+2k$（即$5\mathrm{mod}2=1$）做比较发现十分相似，所以有$x\mathrm{mod}\frac{b}{\gcd (a,b)}=x_1$，此时：
     $$x=\{\begin{array}{l}\frac{b}{\text{gcd}(a,b)},x_1=0\\ (x_1\mathrm{mod}\frac{b}{\text{gcd}(a,b)}+\frac{b}{\text{gcd}(a,b)})\mathrm{mod}\frac{b}{\text{gcd}(a,b)},x_1\ne 0\end{array}$$
     之后根据$y=(c-ax)b$即可求得$y$。

   • 当$y\le 0$时，说明方程没有全是正整数的解，题目要求输出$x_{\min }$和$y_{\min }$，其中$y_{\min }$需要之前求得的最大值$y$做和$x$同样的操作，即补正加判断是否为0，此时即可求得$y_{\min }$。

   • 当$y>0$时，需要输出5个数。此时需要求最小的$x_{\min }$和做过取模处理的$y_{\max }$，以及$x_{\max }$和$y_{\min }$。

     在补正时期就求得$x_{\min }=x$，$y_{\max }=y$，而$y_{\min }$经过$y_{\max }$的取模操作得到，$x_{\max }$通过$\frac{(c-by)}{a}$得到。

     求有多少个正整数解，此时可以通过$\frac{x_{\max }-x_{\min }}{\frac{b}{\gcd (a,b)}}+1$得到（除法得到的是等差数列之间的空隙，需要加1才能得到正整数解）。

#include<bits/stdc++.h>
using namespace std;typedef long long LL;//可改成其他类型，用模板是为了方便拷贝
template<class intT>
intT exgcd(intT a, intT b, intT& x, intT& y) {if (b == 0) {x = 1, y = 0;return a;}intT x1, y1, d;d = exgcd(b, a % b, x1, y1);x = y1, y = x1 - a / b * y1;return d;
}void ac() {LL a, b, c;cin >> a >> b >> c;LL x1, y1, d;d = exgcd(a, b, x1, y1);if (c % d) {cout << -1 << endl;return;}//补正操作//这里需要先算c/d，因为c/d是一个整体，在整数除法顺序变了结果可能会错x1 = c / d * x1;LL g1 = b / d, g2 = a / d;//方便记录b/gcd(a,b)、a/gcd(a,b)x1 = (x1 % g1 + g1) % g1;x1 = x1 == 0 ? g1 : x1;//x_miny1 = (c - a * x1) / b;//y_max//分类讨论if (y1 <= 0) {y1 = (y1 % g2 + g2) % g2;y1 = y1 == 0 ? g2 : y1;cout << x1 << ' ' << y1 << endl;return;}else {//0表示最小值，1表示最大值LL x0 = x1, y0 = y1;//x_miny0 = (y0 % g2 + g2) % g2;y0 = y0 == 0 ? g2 : y0;//y_minx1 = (c - b * y0) / a;//x_maxcout << (x1 - x0) / g1 + 1 << ' ' << x0 << ' ' << y0<< ' ' << x1 << ' ' << y1 << endl;}
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;cin >> T;while (T--)ac();return 0;
}

线性同余方程

同余方程的概念较为复杂，这里只举例算法竞赛常考的线性同余方程$ax\equiv b(\text{mod }m)$。

竞赛中常考求一次同余方程（或线性同余方程） $ax\equiv b(\text{mod }m)$ 的解。

或求方程$ax\mathrm{mod}m=b$，$b\in [0,m-1]$的解。

这个是个人根据实际应用时遇到的问题进行的补充。

例如：求 $4x\equiv 2(\text{mod }6)$，其中一个解为 $x=2$。

【扩欧算法求同余方程的解】

同余方程可以转换为二元一次不定方程：

• 由 $ax\equiv b(\text{mod }m)$ 可得：$ax=ym+b$（$ax÷m=y\dots \dots b$），移项得 $ax-my=b$；

  $ax\equiv b(\text{mod }m)$得到

  $ax÷m=y_1\dots r$，

  $b÷m=y_2\dots r$，$r=b-m{y}_2$，代入$ax÷m=y_1\dots r$有

  $ax=m{y}_1+b-my2=m(y_1-y_2)+b$，为了求$x$，令$y=y_1-y_2$，

  方程变成$ax=ym+b$。

• 最终目的是求$x$，只有少部分情况需要求$y$，所以令 $y=-y$，就变成熟悉的二元一次不定方程：$ax+my=b$；

• 由裴蜀定理得，当 $\text{gcd}(a,m)\mid b$ 时，有解。

那么，就可以用扩展欧几里得算法求出 $x$ 的值。

【扩欧算法求乘法逆元】

• 特殊的，如果 $b$ 等于 1，原方程变为 $ax\equiv 1(\text{mod }m)$，解出的 $x$ 实际上就是 $a$ 在模 $m$ 意义下的乘法逆元。

因此，利用扩欧算法也可以求出乘法逆元。并且不需要模数 $m$ 是质数，仅需 $\text{gcd}(a,m)=1$，也就是 $a,m$ 互质即可。

同余方程OJ示例

扩欧求a的逆元 Nowcoder

【模板】同余方程

真模板题，同余方程$ax\equiv 1(\mathrm{mod}b)$可转换成同余方程$ax+by=1$，当$\text{gcd}(a,b)\ne 1$时方程无解，否则就用扩展欧几里得算法求解。且$x\ne 0$。因为$\text{gcd}(a,b)=1$，$x=0$时，0除以任何数都为0，最大公约数不可能为1。

#include<bits/stdc++.h>
using namespace std;typedef long long LL;LL exgcd(LL a, LL b, LL& x, LL& y) {if (!b) {x = 1, y = 0;return a;}LL x1, y1, d;d = exgcd(b, a % b, x1, y1);x = y1; y = x1 - a / b * y1;return d;
}void ac() {LL a, b;cin >> a >> b;//ax %b = 1->ax+by=1LL x1, y1, d = exgcd(a, b, x1, y1);if (d != 1) {cout << -1 << endl;return;}//d=1，这里可省略x1 = ((x1 / d) % (b / d) + (b / d)) % (b / d);x1 = x1 == 0 ? b : x1;cout << x1 << endl;
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;cin >> T;while (T--)ac();return 0;
}

青蛙的约会 - 洛谷

P1516 青蛙的约会 - 洛谷

1631：【例 1】青蛙的约会

根据题意可得$(x+mt)\equiv (y+nt)(\mathrm{mod}L)$，例如样例：

$(1+3t)\equiv (2+4t)(\mathrm{mod}5)$，当$t=4$时成立。

这里目标是求$t$，所以将方程$(x+mt)\equiv (y+nt)(\mathrm{mod}L)$转化成等式：

$(x+mt)\mathrm{mod}L=(y+nt)\mathrm{mod}L$，因为$x,y\le L$，所以

$x+(mt)\mathrm{mod}L=y+(nt)\mathrm{mod}L$，移项得

$((m-n)t)\mathrm{mod}L=y-x$，

设$a=m-n$，$b=y-x$，则原来的方程变成

$at\mathrm{mod}L=b$（$at\equiv b(\mathrm{mod}L)$），

所以问题变成了求$at+Ly=b$这个方程的最小整数解，$t$的最小整数解就是答案。

但因为过程中出现很多减法，需要先进行补正操作。

#include<bits/stdc++.h>
using namespace std;typedef long long LL;LL exgcd(LL a, LL b, LL& x, LL& y) {if (!b) {x = 1, y = 0;return a;}LL x1, y1, d = exgcd(b, a % b, x1, y1);x = y1; y = x1 - a / b * y1;return d;
}void ac() {LL x, y, m, n, l;cin >> x >> y >> m >> n >> l;//求ax+by=c的解LL a = ((m - n) % l + l) % l, b = l, c = ((y - x) % l + l) % l;LL d = exgcd(a, b, x, y);if (c % d) {cout << "Impossible" << endl;return;}x = x * c / d;x = (x % (b / d) + (b / d)) % (b / d);cout << x << endl;
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;//cin >> T;while (T--)ac();return 0;
}

高次同余方程

问题：给定整数 $a,b,p$，其中 $a,p$ 互质，求一个非负整数 $x$，使得 $a^x\equiv b(\text{mod }p)$。

该问题利用 BSGS(baby step giant step)算法求解。

baby step gaint step，婴儿一小步，巨人一大步。

由欧拉定理的推论知$a^x\equiv a^{x\mathrm{mod}\phi (x)}(\mathrm{mod}p)$，因为$\phi (x)=x\times \frac{p_1-1}{p_1}\times \frac{p_2-1}{p_2}\times \cdots \times \frac{p_k-1}{p_k}$，后面的$\frac{p_i-i}{p_i}$都小于1，所以$\phi (x)<x$，理论上枚举[1,p]就能找到高次同余方程的解。

BSGS算法：

设 $x=im-j$，其中 $m=\lceil \sqrt{p}\rceil ,0\le j\le m-1$，则方程变为 $a^{im-j}\equiv b(\text{mod }p)$。根据同乘性有 $(a^m{)}^i\equiv b\times a^j(\text{mod }p)$。

1. 对于所有的 $j\in [0,m-1]$，把 $b\times a^j\mathrm{mod}p$ 插入一个 Hash 表。

2. 枚举 $i$ 的所有可能取值，即 $i\in [0,m]$，计算出 $(a^m{)}^i\text{ mod }p$，在 Hash 表中查找是否存在对应的 $j$，存在的话更新结果$x=im-j$。其中第1次出现的结果即为最小正整数解。

时间复杂度为 $O(\sqrt{p})$。

需要处理很多边界条件，具体在代码中实现。至于为什么会有这些边界条件，等本人到达那个境界再补充。

下面给出利用 BSGS 求解同余方程 $a^x\equiv b\text{ mod }p$ 的最小非负整数解，无解时输出 $-1$。

【代码实现】

typedef long long LL;
//求a^x % p == b的解
LL BSGS(LL a, LL b, LL p) {a %= p; b %= p;if (a == 0)//a%p = b%p=0, 方程有无限多个解return b == 0 ? 1 : -1;if (b == 1)//a^0 %1 =1成立return 0;//BSGS要求a、p互质，否则a模p的乘法逆元可能不存在//扩展为Pohlig-Hellman 算法 或先处理公共因子LL x, y, d = exgcd(a, p, x, y);if (d > 1) {if (b % d)return -1;//方程转化为 (a/d)^x ≡ (b/d) * (a/d)^{-1} (mod p/d)p /= d;d = exgcd(a / d, p, x, y);if (d != 1)return -1;x = (x % p + p) % p;x = x == 0 ? p : x;//求a/d的逆元b = qmul(b / d, x, p);return BSGS(a, b, p) + 1;}map<LL, LL>mp;//键值对LL m = sqrt(p) + 1;//预处理b^jfor (LL j = 0, step = b; j < m; j++) {if (!mp.count(step))mp[step] = j;step = qmul(step, a, p);}//预处理a^mLL am = 1;for (LL i = 1; i <= m; i++)am = qmul(am, a, p);for (LL i = 1, step = 1; i <= m; i++) {step = qmul(step, am, p);if (mp.count(step))return i * m - mp[step];}return -1;
}

高次同余方程OJ示例

模板题计算器 - 洛谷

[P2485 SDOI2011] 计算器 - 洛谷

模板题，但高次同余方程需要规避很多边界条件。

#include<bits/stdc++.h>
using namespace std;typedef long long LL;LL exgcd(LL a, LL b, LL& x, LL& y) {if (!b) {x = 1; y = 0;return a;}LL x1, y1, d = exgcd(b, a % b, x1, y1);x = y1; y = x1 - a / b * y1;return d;
}LL qmul(LL a, LL b, LL m) {LL ans = 0;while (b) {if (b % 2)ans = (ans + a) % m;a = (a * 2) % m;b /= 2;}return ans;
}//求a^x % p == b的解
LL BSGS(LL a, LL b, LL p) {a %= p; b %= p;if (a == 0)//a%p = b%p=0, 方程有无限多个解return b == 0 ? 1 : -1;if (b == 1)//a^0 %1 =1成立return 0;//BSGS要求a、p互质，否则a模p的乘法逆元可能不存在//扩展为Pohlig-Hellman 算法 或先处理公共因子LL x, y, d = exgcd(a, p, x, y);if (d > 1) {if (b % d)return -1;//方程转化为 (a/d)^x ≡ (b/d) * (a/d)^{-1} (mod p/d)p /= d;d = exgcd(a / d, p, x, y);if (d != 1)return -1;x = (x % p + p) % p;x = x == 0 ? p : x;//求a/d的逆元b = qmul(b / d, x, p);return BSGS(a, b, p) + 1;}map<LL, LL>mp;//键值对LL m = sqrt(p) + 1;//预处理b^jfor (LL j = 0, step = b; j < m; j++) {if (!mp.count(step))mp[step] = j;step = qmul(step, a, p);}//预处理a^mLL am = 1;for (LL i = 1; i <= m; i++)am = qmul(am, a, p);for (LL i = 1, step = 1; i <= m; i++) {step = qmul(step, am, p);if (mp.count(step))return i * m - mp[step];}return -1;
}void ac(LL k) {LL y, z, p;cin >> y >> z >> p;if (k == 1) {//快速幂LL ans = 1;while (z) {if (z % 2)ans = qmul(ans, y, p);y = qmul(y, y, p);z /= 2;}cout << ans << endl;}else if (k == 2) {//求同余方程y*x+p*t=zLL x, t, d = exgcd(y, p, x, t);if (z % d) {cout << "Orz, I cannot find x!\n";return;}x = qmul(x, z / d, p / d);x = (x % p + p) % p;cout << x << endl;}else {LL ans = BSGS(y, z, p);if (ans == -1)cout << "Orz, I cannot find x!\n";elsecout << ans << endl;}
}int main() {ios::sync_with_stdio(false);cin.tie(0);LL T, K;cin >> T >> K;while (T--)ac(K);return 0;
}

模板题BSGS - 洛谷

[P3846 TJOI2007] 可爱的质数/【模板】BSGS - 洛谷

还是模板题。

#include<bits/stdc++.h>
using namespace std;typedef long long LL;LL exgcd(LL a, LL b, LL& x, LL& y) {if (!b) {x = 1; y = 0;return a;}LL x1, y1, d = exgcd(b, a % b, x1, y1);x = y1; y = x1 - a / b * y1;return d;
}LL qmul(LL a, LL b, LL m) {LL ans = 0;while (b) {if (b % 2)ans = (ans + a) % m;a = (a * 2) % m;b /= 2;}return ans;
}//求a^x % p == b的解
LL BSGS(LL a, LL b, LL p) {a %= p; b %= p;if (a == 0)//a%p = b%p=0, 方程有无限多个解return b == 0 ? 1 : -1;if (b == 1)//a^0 %1 =1成立return 0;//BSGS要求a、p互质，否则a模p的乘法逆元可能不存在//扩展为Pohlig-Hellman 算法 或先处理公共因子LL x, y, d = exgcd(a, p, x, y);if (d > 1) {if (b % d)return -1;//方程转化为 (a/d)^x ≡ (b/d) * (a/d)^{-1} (mod p/d)p /= d;d = exgcd(a / d, p, x, y);if (d != 1)return -1;x = (x % p + p) % p;x = x == 0 ? p : x;//求a/d的逆元b = qmul(b / d, x, p);return BSGS(a, b, p) + 1;}map<LL, LL>mp;//键值对LL m = sqrt(p) + 1;//预处理b^jfor (LL j = 0, step = b; j < m; j++) {if (!mp.count(step))mp[step] = j;step = qmul(step, a, p);}//预处理a^mLL am = 1;for (LL i = 1; i <= m; i++)am = qmul(am, a, p);for (LL i = 1, step = 1; i <= m; i++) {step = qmul(step, am, p);if (mp.count(step))return i * m - mp[step];}return -1;
}int main() {ios::sync_with_stdio(false);cin.tie(0);LL p, b, n;cin >> p >> b >> n;p = BSGS(b, n, p);if (p == -1)cout << "no solution";elsecout << p << endl;return 0;
}

再看乘法逆元

【乘法逆元】

乘法逆元，一般用于求 $\frac{a}{b}(\mathrm{mod}p)$ 的值。

先算出 $b$ 在模 $p$ 意义下的乘法逆元 $b^{-1}$，然后计算 $a\times b^{-1}(\text{mod }p)$ 的值，就可以将除法转化为乘法。即$\frac{a}{b}\equiv a\times b^{-1}(\mathrm{mod}p)$。

求单个元素模m的逆元

【方法一：费马小定理 + 快速幂】

• 问题：求 $a$ 在模 $m$ 意义下的逆元。

• 前置知识回顾：

  • 费马小定理：如果 $p$ 为质数，且 $a,p$ 互质，则 $a^{p-1}\equiv 1(\mathrm{mod}p)$。
    且$a\times a^{p-2}\equiv 1(\mathrm{mod}p)$，则$a^{p-2}$就是$a$的模$p$乘法逆元。
  • 快速幂：利用二进制和模运算特性求$a^b\mathrm{mod}m$的算法。

• 前提条件：$a,m$ 互质，且 $m$ 是质数；

• 方法：求出 $a^{m-2}\mathrm{mod}m$，就是 $a$ 在模 $m$ 意义下的逆元。

• 时间复杂度：$O(\log m)$。

【方法二：扩欧算法求逆元】

• 问题：求 $a$ 在模 $m$ 意义下的逆元。

• 前置知识回顾：

  • 扩展欧几里得算法的结论：对不定方程$ax+by=\text{gcd}(a,b)$，可通过一下结论和递归求得不定方程的特解：

    $\{\begin{array}{l}x=y_1,\\ y=x_1-\frac{a}{b}{y}_1\end{array}$

    其中递归的边界条件是$b=0$，此时$x_1=1,y_1=0$。

• 前提条件：$a,m$ 互质；

• 方法：解同余方程 $ax\equiv 1(\mathrm{mod}m)$（或$ax\mathrm{mod}m=1$）。将方程转换成$ax+my=1$，求方程的解中$x$可能得值。

• 时间复杂度：$O(\log m)$。

递推求前n个数模p的逆元

• 问题：求 $1\sim n$ 中所有的数在模 $p$ 意义下的逆元。

• 前提条件：$1\sim n$ 中所有的数都与 $p$ 互质。最常见的情形是$p>n$且$\text{gcd}(n,p)=1$。

• 方法：线性递推：

  对$1\sim n\mathrm{mod}p$：

  • 当 $i=1$：inv[1] = 1。

  • 当 $i>1$：inv[i] = p - (p/i) * inv[p % i] % p。

    设 i ∈ { 1 , 2 , … , n } i\in \{1,2,\dots,n\} i∈{1,2,…,n}，则$p÷i=k\dots r$，

    于是$p=i\times k+r$，

    所以有$(i\times k+r)\equiv 0(\mathrm{mod}p)$，或$(i\times k+r)\mathrm{mod}p=0$。

    两边同时乘$i^{-1}\mathrm{mod}p$和$r^{-1}\mathrm{mod}p$有

    $(k{r}^{-1}+i^{-1})\equiv 0(\mathrm{mod}p)$（$(k{r}^{-1}+i^{-1})\mathrm{mod}p=0$），

    所以$i^{-1}\mathrm{mod}p=-k{r}^{-1}\mathrm{mod}p$，即$i^{-1}\equiv (-k{r}^{-1})(\mathrm{mod}p)$

    $k$和$r^{-1}\mathrm{mod}p$用$i$和$p$代替有

    $i^{-1}\equiv (-\frac{p}{i}\times (p\mathrm{mod}i{)}^{-1})(\mathrm{mod}p)$，或

    $i^{-1}\mathrm{mod}p=(-\frac{p}{i}\times (p\mathrm{mod}i{)}^{-1})\mathrm{mod}p$，为防止出现负数，进行补正操作：

    $i^{-1}\mathrm{mod}p=p-((\frac{p}{i}\times (p\mathrm{mod}i{)}^{-1})\mathrm{mod}p)$

    所以根据这个递推式即可求得$1\sim n$ 中所有的数在模 $p$ 意义下的逆元。

• 时间复杂度：$O(n)$。

乘法逆元OJ示例

递推求前n个数模p的逆元

P3811 【模板】模意义下的乘法逆元 - 洛谷

套用递推式即可。但输入、输出的数据量比较大，需要对cin进行优化。

#include<bits/stdc++.h>
using namespace std;typedef long long LL;void ac() {LL n, p;cin >> n >> p;vector<LL>inv(n + 1, 0);inv[1] = 1;for (LL i = 2; i <= n; i++) {//2个都可以inv[i] = p - ((p / i * inv[p % i]) % p);//inv[i] = p - p / i * inv[p % i] % p;//优先级：(* /) % (+ -)}for (LL i = 1; i <= n; i++)cout << inv[i] << '\n';
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;//cin >> T;while (T--)ac();return 0;
}

求单个元素模p的逆元

【模板】逆元

题目没保证$p$是质数，因此不能用费马小定理。所以用扩展欧几里得定理求解。

这题和【模板】同余方程可以共用同一个AC的代码，准确来说是同一个题。

#include<bits/stdc++.h>
using namespace std;typedef long long LL;//扩展欧几里得算法求不定方程的特解
LL exgcd(LL a,LL b,LL&x, LL&y){if(!b){x=1,y=0;return a;}LL x1,y1,d=exgcd(b,a%b,x1,y1);x=y1,y=x1-a/b*y1;return d;
}void ac() {LL a,p;cin>>a>>p;//求逆元转化为求ax+py=1的解LL x,y,d=exgcd(a,p,x,y);if(d!=1){cout<<-1<<'\n';return;}x=(x+p)%p;cout<<x<<'\n';
}int main() {ios::sync_with_stdio(false);cin.tie(0);int T = 1;cin >> T;while (T--)ac();return 0;
}