话说好久没来博客上面写过东西了，之前集训过于辛苦了，但有很大的收获，我觉得有必要把它们拿出来总结分享。之前一直是个数论渣（小学初中没好好念过竞赛的缘故吧），经过一道题目对一些基础算法有了比较深刻的理解，在这里我打算系统地讲出这道题目涉及的大部分内容，希望可以帮到大家。
原题地址：http://acm.sgu.ru/problem.php?contest=0&problem=261
题目大意：给出质数$p$、$k$和一个自然数$a$，求关于$x$的同余方程$x^k \equiv a \pmod p$在区间$[0,\,p-1]$内所有解
数据范围：$2 \le p \le 10^9, 2 \le k \le 100000$
这种同余方程叫做“k次剩余”，这题也是一个模板题，但它涉及了比较多的数论算法，我会一一在这里尽可能讲清楚，这些概念和算法包括：
原根、缩系、离散对数、乘法逆元、线性同余方程、快速幂算法、扩展欧几里德算法、Baby Step - Gaint Step算法
虽然我希望近可能地系统化讲解，但是以下基本概念不会做任何讲解，不懂的同学可以查到大量的相关资料，它们包括：
素数、带余除法、模运算、幂运算、欧拉函数$\varphi(n)$、渐进复杂度分析
首先我们通过引入一些概念对所求进行一些变形，最后再讲述求解所需的算法，前面的理论定义部分可能较为枯燥：
阶、原根、缩系与离散对数：
定义 1（阶）：若$gcd(a,\,m)=1$，满足$a^r \equiv 1 \pmod m$的最小正整数$r$称为$a$模$m$的阶
性质 1（阶的判定定理）：$gcd(a,\,m)=1$时，$a$模$m$的阶是$r$等价于以下两点：
1：$a^r \equiv 1 \pmod m$
2：对于$r$的每一个质因子$p$，均有$a^{\frac{r}{p}} \not\equiv 1 \pmod m$
证明：第一条是阶的定义，第二条应用反证法，若存在质因子$p$使得其不成立，则与$r$的最小性矛盾
定义 2（原根）：若整数$g$模$m$的阶为$\varphi(m)$，则称$g$为$m$的原根
性质 2（原根的存在情况）：正整数$m$存在原根，当且仅当$m=2,\,4,\,p^a,\,2p^a$，其中$p$是奇素数且$a \ge 1$
由性质2我们可以得知所有的质数都有原根。
性质 3（原根的分布）：一个数可能对应多个原根，而且最小的原根一般都很小（基本没例外）
定义 3（缩系）：若模$m$有原根$g$，则$\{g^0,\,g^1,\,\dots,\,g^{\varphi(m)-1}\}$称为模$m$的缩系。显然对于每个与$m$互质的正整数$a$，在$m$的缩系中仅存在唯一的一个整数$k$，使得$g^r \equiv a \pmod m$
定义 4（离散对数）：我们将形如$a^k \equiv n \pmod m$（m是质数）的式子中的$k$称为$n$以$a$为底在模$m$下的离散对数，特别地，在$a$为$m$的原根$g$时，我们称$k$为$n$模$m$的离散对数，由于$k$的值与原根$g$有关，我们将它记作$ind_{g}n$
由于我们讨论的模是质数，所以我们可以将不超过$m$的运算全部放到它的缩系里面来做了。而且不超过$m$的整数$n$模$m$的离散对数正巧与$m$的缩系中$n$对应的值相等。
性质 4（离散对数运算性质）：
类似于对数的运算性质，离散对数的性质大致有以下几条：
1：若$gcd(a,\,m)=gcd(b,\,m)$，$a \equiv b \pmod m$等价于$ind_{g}a=ind_{g}b$
2：$ind_{g}ab=ind_{g}a+ind_{g}b \pmod{\varphi(m)}$
3：$ind_{g}a^n=n\,ind_{g}a \pmod{\varphi(m)}$
以上涉及的内容全在整数范围内
好了，进行了那么多定义和基础知识的灌输，题目涉及的理论基础也差不多说完了，让我们来回归题目本身吧：
$x^k \equiv a \pmod p$等价于$k\,ind_{g}x \equiv ind_{g}a \pmod{\varphi(p)}$，其中$g$是模$p$的原根（这个过程强烈建议读者手推）。这样，我们直接求解这个线性同余方程就好了，问题被分解为了求原根、离散对数以及求解线性同余方程，接下来我们开始讲解每步的求解过程，可以选择不会的地方跳着阅读
求解原根：
由上面的性质3可以知道最小的原根一般都很小，由于目前不存在寻找原根的比暴力枚举更有效的方法，所以我们直接枚举原根然后进行判定就好了
算法流程：
1.输入一个质数$p$
2.枚举1到$p-1$的所有整数$k$，判断$p-1$是否为$k$的阶（即$k^{p-1}$是否同余于1，并枚举$p-1$的所有质因子$r$判断$k^{\frac{p-1}{r}}$是否不同余于1，若有任何一条不满足，则$k$不是原根，如果两条都满足，则返回$k$）
算法代码：
inline long long findRoot(long long p)
{
long long s = ceil(sqrt(p - )); for(long long i = ; i < p; ++i)
{
if(pow(i, p - , p) == 1long long)
{
long long sign = ;
for(long long j = ; j <= s; ++j)
if(((p - ) % j == ) && (pow(i, (p - ) / j, p) == ))
{
sign = ; break;
}
if(sign)return i;
}
}
}
快速幂：
上面寻找原根的程序中调用了一个函数pow(n, k, p)，用来计算$n^k\mod p$的值。这个算法用到了指数运算的一些简单技巧
算法流程：
1.若k = 0, 返回 1
2.设w = pow(n*n mod p, k / 2, p)
3.若k为奇数，返回w * n mod p, 否则返回w
这个算法过于简单了，不多说
时间复杂度：$O(\lg k)$
算法代码：
long long pow(long long n, long long k, long long p)
{
if(k == 0long long)return 1long long;
long long w = ;
if(k & )w = n % p;
long long ans = pow(n * n % p, k >> , p);
return ans * w % p;
}
接下来的工作就是求离散对数啦：我们采用Baby Step-Gaint Step算法（小步大步）
Baby Step - Gaint Step算法：
这一算法在$O(\sqrt{n})$的时间内解决求解离散对数的问题（关于$k$的方程$g^k \equiv n \pmod p$的解）
首先我们设$s=\sqrt{\lceil p \rceil}$，由于$k$的范围始终是$[0,\,p-1]$，所以可以设$k=st+r$
预先计算出所有小于等于$s$的$i$的$g^i$的值，然后从0到$s$枚举$t$，之后检查是否有满足条件的$r < s$存在，如果有，则$ts+r$则是符合条件的解。
注意寻找是否有满足条件的$r$的时候，我们有两种处理方式：第一种是使用哈希表或树形结构维护，第二种是将它们预处理之后排序再二分。这里为了方便我使用了标准库中的map
算法代码：
inline long long ind(long long x, long long p, long long g)
{
static map<long long, long long> tmp;
long long s = ceil(sqrt(p));
long long w = pow(g, s, p); for(long long r = ; r <= s; ++r)tmp.insert(make_pair(pow(g, r, p), r));
for(long long t = ; t <= s; ++t)
{
long long gt = pow(w, t, p);
long long anti = findv(gt);
if(tmp.count(x * anti % p))return (tmp[x * anti % p] + t * s);
}
return -;
}
注意到在这个算法中，查找合条件的$r$时我们使用了除法，但是在模意义下除法会出问题，我们该为乘以除数模意义下的乘法逆元。
乘法逆元：
定义 5（乘法逆元）：若$xk \equiv 1 \pmod p$，则称$k$和$x$互为模$p$意义下的乘法逆元。
解法：将同余方程$xk \equiv 1 \pmod p$转化为等式方程$xk-np=1\quad(n \in N)$，使用下面的扩展欧几里德算法求解
扩展欧几里德算法：
这个是比较基础但是讲起来比较麻烦的数论算法了，这里只讲推导过程和做法，严谨的证明请查阅相关书籍。
此算法用于求解形如$ax+by=gcd(a,\,b)$的不定方程的解$(x,\,y)$。根据欧几里德算法我们知道有$bx_0+(a \mod b)y_0=gcd(a,\,b)$
经过一些整理我们可以发现由于$a \mod b = a - b \times \lfloor \frac{a}{b} \rfloor$，$y_0a+\left(x_0+y_0 \times \lfloor \frac{a}{b} \rfloor \right)b=gcd(a,\,b)$
所以有$x=y_0, y = \left(x_0+y_0 \times \lfloor \frac{a}{b} \rfloor \right)$，按照欧几里德算法的方式递归下去就好了，这样我们解决了求解逆元的问题。而这道题的最后一步求解同余方程也要使用扩展欧几里德算法解决。
算法代码：
long long exgcd(long long a, long long b, long long &x1, long long &y1)
{
if(b == ){x1 = ; y1 = ; return a;}
long long x0, y0;
long long G = exgcd(b, a % b, x0, y0);
x1 = y0; y1 = x0 - a / b * y0;
return G;
}
求解逆元的代码如下：
inline long long findv(long long x)
{
long long d, t;
exgcd(x, p, d, t);
while(d < )d += p;
return d;
}
（今晚太累了所以后面讲的可能不是很清楚，同学们有不明白的地方可以查阅相关资料或者给我留言，我会尽力完善）
到这里我们总结了一些常用的模运算下的算法，并完整地讲完了这道题，最后附上这道题的完整代码以供查错等需要
//date 20140323
#include <cstdio>
#include <cstring>
#include <cmath>
#include <map>
#include <algorithm>
#include <iostream> using namespace std; typedef long long LL; LL x, k, a, p; LL pow(LL n, LL k, LL p)
{
if(k == 0LL)return 1LL;
LL w = ;
if(k & )w = n % p;
LL ans = pow(n * n % p, k >> , p);
return ans * w % p;
} inline LL findRoot(LL p)
{
LL s = ceil(sqrt(p - )); for(LL i = ; i < p; ++i)
{
if(pow(i, p - , p) == 1LL)
{
LL sign = ;
for(LL j = ; j <= s; ++j)
if(((p - ) % j == ) && (pow(i, (p - ) / j, p) == ))
{
sign = ; break;
}
if(sign)return i;
}
}
} LL exgcd(LL a, LL b, LL &x1, LL &y1)
{
if(b == ){x1 = ; y1 = ; return a;}
LL x0, y0;
LL G = exgcd(b, a % b, x0, y0);
x1 = y0; y1 = x0 - a / b * y0;
return G;
} inline LL findv(LL x)
{
LL d, t;
exgcd(x, p, d, t);
while(d < )d += p;
return d;
} inline LL ind(LL x, LL p, LL g)
{
static map<LL, LL> tmp;
LL s = ceil(sqrt(p));
LL w = pow(g, s, p); for(LL r = ; r <= s; ++r)tmp.insert(make_pair(pow(g, r, p), r));
for(LL t = ; t <= s; ++t)
{
LL gt = pow(w, t, p);
LL anti = findv(gt);
if(tmp.count(x * anti % p))return (tmp[x * anti % p] + t * s);
}
return -;
} LL ans[], rans[];
int nans; int main()
{
freopen("sgu.in", "r", stdin);
freopen("sgu.out", "w", stdout); cin >> p >> k >> a; if(a >= p)
{
printf("0\n");
return ;
} if(a == )
{
printf("1\n0\n");return ;
}
LL g = findRoot(p);
LL q1, q2 = ind(a, p, g), w; if(w == -)
{
printf("0\n");
return ;
} LL G = exgcd(k, p - , q1, w); if(q2 % G != )
{
printf("0\n");
return ;
} q1 = q1 * (q2 / G) % (p - );
LL delta = (p - ) / G; nans = ;
for(int i = ; i < G; ++i)
{
q1 = ((q1 + delta) % (p - ) + p - ) % (p - );
ans[++nans] = pow(g, q1, p);
} sort(ans + , ans + nans + ); int tot = ;
for(int i = ; i <= nans; ++i)
{
if(ans[i] != ans[i - ])rans[++tot] = ans[i];
}
cout << tot << endl;
for(int i = ; i < tot; ++i)cout << rans[i] << ' ';
cout << rans[tot] << endl;
return ;
}
SGU 261
结束语：数论在计算机科学中占据着重要的地位，这道简单的题也只能算是管中窥豹，后面我会学习新的知识并于大家分享