前言
难以置信, 遇到这个问题竟然是因为做一道CTF Reverse的题, 实际上我以前打ACM的时候都没有写过任意模数的K次剩余这种东西. 最简单的模质数的K次剩余的例题是
“`HDU 3930 Broot“`, 进阶版本(模任意奇数的K次剩余)的例题是“`BZOJ 2219 数论之神“`, 而终极版本(任意模数的K次剩余)似乎可以交的题只有“`51nod 1123“`, 截至本文写作时只有44个AC. 然而这个CTF题本身就有50多个通过了. ~~CTF>ACM 实锤~~ (暴论)
其实是因为这个Reverse的题里面幂指数是质数，所以会跟模数互质，然后就可以用RSA的方法来做了. 其次该题目更多的工作在于去花指令混淆上, 而非算法破解.
OK那么就来看一下这个问题从最简单的版本到终极版本怎么求解吧…
模质数的K次剩余
已知$a$, $b$, $p$, 求使得
$$
x^{a} \equiv b \quad(\bmod \ p)
$$
成立的所有 $x$. 其中 $p$ 是质数.
由于 $p$ 是质数, 所以 $p$ 存在原根 $g$ , 此时对于模 $p$ 意义下的任意数 $w$ ($0\le w \le p-1$) 存在唯一的 $i$ ($0\le i \le p-1$) 使得 $w\equiv g^i\quad(\bmod \ p)$.
由此可以将最终的答案用 $g$ 来表示: $x=g^{c}$, 带入上式转化为求解
$$
(g^{c})^{a}\equiv (g^{a})^c\equiv b \quad(\bmod \ p)
$$
此时 $g$, $a$, $b$, $p$, 已知, 只需要解出来 $c$. 此时相当于求解离散对数, 使用 Baby-Step-Giant-Step 可以在 $\mathcal{O}(\sqrt{p})$ 时间内得到一个特解 $x_0 \equiv g^{c} \quad(\bmod \ p)$.
在已知一个特解的情况下求出所有解是简单的, 由原根的性质可知 $g^{\varphi(p)}\equiv 1 \quad(\bmod \ p)$, 因此:
$$
\forall t \in \mathbb{Z}, x^{a} \equiv g^{c \cdot a+t \cdot \varphi(p)} \equiv b \quad(\bmod \ p)
$$
因此所有的解为:
$$
\forall t \in \mathbb{Z}, a \mid t \cdot \varphi(p), x \equiv g^{c+\frac{t \cdot \varphi(p)}{a}} \quad(\bmod p)
$$
上面幂次部分要能整除必须要有$\frac{a}{\operatorname{gcd}(a, \varphi(p))} \mid t$, 可以设 $t=\frac{a}{\operatorname{gcd}(a, \varphi(p))} \cdot i$, 于是得到所有的解为:
$$
\forall i \in \mathbb{Z}, x \equiv g^{c+\frac{\varphi(p)}{\operatorname{gcd}(a, \varphi(p))} \cdot i} \quad(\bmod p)
$$
HDU 3930 Broot 代码:
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef __int128 i128;
const int maxn=100;
/*
* Complexity:
* Get primitive root: O(sqrt(M)+0.25*M*log(M))
* Discrete logarithm: O(sqrt(M))
* Total: O(sqrt(M)+0.25*M*log(M))
*/
ll gcd(ll a, ll b){
return b?gcd(b,a%b):a;
}
ll qpow(i128 a, i128 b, i128 p){
ll res=1;
while(b){
if(b&1)
res=(res*a)%p;
a=(a*a)%p;
b>>=1;
}
return res;
}
int pi[maxn],pcnt;
ll get_g(ll P){
pcnt=0;
ll phi= P - 1;
ll tphi=phi;
if(!(tphi & 1ll)){
pi[++pcnt]=2;
while(!(tphi & 1)) tphi>>=1;
}
for(ll i=3;i*i <= tphi; i+=2){
if(tphi % i == 0){
pi[++pcnt]=i;
while(tphi % i == 0) tphi/=i;
}
}
if(tphi != 1) pi[++pcnt]=tphi;
for(ll g=2; g < P; ++g){
bool flag=true;
for(int i=1;i<=pcnt;++i){
if(qpow(g,phi/pi[i], P) == 1){
flag=false;
break;
}
}
if(flag) return g;
}
return -1;
}
unordered_map<ll,int> Hash;
ll bsgs(ll a, ll b, ll M){
if(b==1) return 0;
Hash.clear();
ll sq=sqrt(M*1.0);
i128 x=1,p=1;
for(int i=0;i<sq;++i,p=p*a%M) Hash[p*b%M]=i;
for(ll i=sq,j;i<=M;i+=sq){
if(Hash.find(x=x*p%M)!=Hash.end()){
return i-Hash[x];
}
}
return -1;
}
vector<ll> ans;
void k_residue(ll k, ll y, ll M){
ans.clear();
if(y==0){
ans.push_back(0);return ;
}
ll g=get_g(M);
ll c=bsgs(qpow(g,k,M),y,M);
if(c==-1) {
ans.push_back(-1);return;
}
ll delta=(M-1)/gcd(k,M-1);
for(ll cur=c%delta;cur<M-1;cur+=delta)
ans.push_back(qpow(g,cur,M));
return ;
}
int main() {
ll k,y,M;
int cs=0;
while(scanf("%lld %lld %lld",&k,&M,&y)==3){
++cs;printf("case%d:\n",cs);
k_residue(k,y,M);
sort(ans.begin(),ans.end());
for(int i=0,len=ans.size();i<len;++i){
printf("%lld\n",ans[i]);
}
}
return 0;
}
模任意奇数的K次剩余
已知$a$, $b$, $p$, 求使得
$$
x^{a} \equiv b \quad(\bmod \ p)
$$
成立的所有 $x$. 其中 $p$ 是任意奇数.
此时套用原来的解法会存在以下两个问题:
- $p$ 不一定有原根
- 即使 $p$ 有原根, $b$ 跟 $p$ 不互质的话则无法用原根表示 $b$, 也就无法转化为离散对数求解.
对于问题1, 可以考虑将 $p$ 做质因子分解, 即$p=p_{1}^{e_1}…p_{j}^{e_j}$. 由原根的存在性定理可知 $p_{j}^{e_j}$ 一定有原根. 将每个 $p_{j}^{e_j}$ 作为模数分别求解得到 $n_j$ 个解, 求完了之后对于其中的每一个解, 令其与模数为 $p_{1}^{e_1}…p_{j-1}^{e_{j-1}}$ 时的所有解(假设有 $n_{j-1}$ 个)用中国剩余定理两两合并, 这可以用两个数组迭代地实现. 本次合并完后应该得到 $n_{j}\times n_{j-1}$ 个解, 注意每次合并完之后去重一下.
对于问题2, 我们在将问题分解到对每个 $p_{j}^{e_j}$ 上之后, 对 $b \ \bmod(\ p_{j}^{e_j})$ 的结果与 $p_{j}^{e_j}$ 的关系进行分类讨论.
- $gcd(b \ \bmod(\ p_{j}^{e_j}), p_{j}^{e_j})=1$
此时 $b \ \bmod(\ p_{j}^{e_j})$ 可以表示为 $g^c$, 用数论的语言来说就是可以求指标. 此时直接套用模质数的K次剩余即可.
-
$1\lt gcd(b \ \bmod(\ p_{j}^{e_j}), p_{j}^{e_j}) \lt p_{j}^{e_j}$
假设 $x^a \ \bmod(\ p_{j}^{e_j})=b \ \bmod(\ p_{j}^{e_j})=u\times p_{j}^{k}$, 其中 $gcd(u,p)=1$. 令同余等式同时除以一个 $p_{j}^{k}$ 得到 $(\frac{x}{p_{j}^{\frac{k}{a}}})^a \ \bmod(\ p_{j}^{e_j-k})=u$ (如果不能除说明无解). 此时模数$p_{j}^{e_j-k}$ 与 $u$ 互质, 回到情况1, 直接套用模质数的K次剩余, 在求出来特解之后, 找所有解的过程中将每个解乘上 $p_{j}^{\frac{k}{a}}$得到原问题的解, 并且乘上之后要对$p_{j}^{k}$ 取模(之前的套用模质数的时候运算都需要对$p_{j}^{e_j-k}$取模, 包括枚举解过程中的增量).
-
$gcd(b \ \bmod(\ p_{j}^{e_j}), p_{j}^{e_j}) = p_{j}^{e_j}$
此时其实意味着$b \ \bmod(\ p_{j}^{e_j})=0$. 若$x^a\equiv 0 \ \bmod(\ p_{j}^{e_j})$, 则必有 $p^{\lceil\frac{e_j}{a}\rceil}|x$. 此时只需要枚举其倍数即可, 则解为 $\forall i \in \mathbb{Z}, x \equiv i\times p^{\lceil\frac{e_j}{a}\rceil} \ \bmod(\ p_{j}^{e_j})$
这样就解决了模任意奇数K次剩余的问题, 例题是
“`BZOJ 2219 数论之神“`, 我没有写这个题的代码(BZOJ已经无了), 而是直接写了更进阶的任意模数K次剩余的代码, 见下一小节.
模任意数的K次剩余
已知$a$, $b$, $p$, 求使得
$$
x^{a} \equiv b \quad(\bmod \ p)
$$
成立的所有 $x$. 其中 $p$ 是任意整数.
此时再套用模奇数的K次剩余的方法, 有一个问题是 $p$ 的因子中可能会包括 $2$ , 而我们知道对于 $\forall k \in \mathbb{Z},2^k$ , 只有 $2$ 和 $4$ 是有原根的:
原根存在定理：一个数 $m$ 存在原根当且仅当 $m=2,4, p^{\alpha}, 2 p^{\alpha}$, 其中 $p$ 为奇素数, $\alpha \in \mathbb{N}^{*}$
那么我们在将模数分割到每个 $p_j^{e_j}$ 上的时候就需要对 $p_j=2$ 这种情况做特殊处理. 此时不能用原根将问题转化为离散对数问题了, 前面的分析全部失效, 但当解的数量较小的时候, 可以用倍增的方法来求出所有模 $2^k$ 意义下的解.
具体地, 假设我们要求的最终答案为 $x$, 满足 $x^{a} \equiv b \quad(\bmod \ p)$, 且 $gcd(2^e,p)=2^e, e>0$, 那么如果上式在模 $2^{k-1}$ 意义下存在一个解为 $x_i$, 则它对应于 $2^k$ 的解只能是 $x_i$ 或者 $x_{i}+2^{k-1}$. 我们只需要对于每个$i$把这两个可能的解代进$x^{a} \equiv b \quad(\bmod \ 2^i)$里检验一下就可以了, 如果不能使得式子成立就丢掉. 当然直观来看这种做法似乎会得到最多 $2^{e}$ 个解, 但实际上一般解是很少的, 可以保证较快的算出来.
将模 $2^e$ 意义下的所有解求出之后, 其他的就和模任意质数的K次剩余一样了, 对不同的质因子求一遍然后用中国剩余定理合并即可. 这样的话就解决了任意模数K次剩余的问题.
时间复杂度:
- 质因子分解 $\mathcal{O}(\sqrt{p})$
- 找原根 $\mathcal{O}(\sqrt{p}+0.25\times p\times log(p))$
- 离散对数BSGS $\mathcal{O}(\sqrt{p})$
- CRT合并, $2^e$ 因子的处理 $\mathcal{O}(答案的个数*log(p))$
例题
“`51nod 1123 X^A mod B“`, 代码:
#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef __int128 i128;
const int maxn=100;
ll upperInt(ll a, ll b){
return (a+b-1)/b;
}
ll gcd(ll a, ll b){
return b?gcd(b,a%b):a;
}
ll extEuclid(ll a,ll b,ll &x,ll &y){
ll d;
if(b==0){
x=1;y=0;
return a;
}
d=extEuclid(b,a%b,y,x);
y-=a/b*x;
return d;
}
ll qpow(i128 a, i128 b, i128 p){
ll res=1;
while(b){
if(b&1)
res=(res*a)%p;
a=(a*a)%p;
b>>=1;
}
return res;
}
void factor(ll P, ll pi[], ll ai[], ll &pcnt){
pcnt=0;
ll tP=P;
if(!(tP & 1ll)){
pi[++pcnt]=2;
ai[pcnt]=0;
while(!(tP & 1)) tP>>=1,ai[pcnt]++;
}
for(ll i=3;i*i <= tP; i+=2){
if(tP % i == 0){
pi[++pcnt]=i;
ai[pcnt]=0;
while(tP % i == 0) tP/=i,ai[pcnt]++;
}
}
if(tP != 1) {
pi[++pcnt]=tP;
ai[pcnt]=1;//error2
}
}
ll piG[55],aiG[55],pcntG;
ll getG(ll P, ll phi){
factor(phi,piG,aiG,pcntG);
for(ll g=2; g < P; ++g){
bool flag=true;
for(int i=1;i<=pcntG;++i){
if(qpow(g,phi/piG[i], P) == 1){
flag=false;
break;
}
}
if(flag) return g;
}
return -1;
}
unordered_map<ll,int> Hash;
ll bsgs(ll a, ll b, ll M){
if(b==1) return 0;
Hash.clear();
ll sq=sqrt(M*1.0);
i128 x=1,p=1;
for(int i=0;i<sq;++i,p=p*a%M) Hash[p*b%M]=i;
for(ll i=sq;;i+=sq){
if(Hash.find(x=x*p%M)!=Hash.end()){
return i-Hash[x];
}
if(i>M) break;//error3
}
return -1;
}
ll as[55],ps[55];
ll crt(ll as[], ll ps[], ll len)
{
int i;
long long d,x,y,m,n,ret;
ret=0;
n=1;
for (i=0;i<len;i++)
n*=ps[i];
for (i=0;i<len;i++){
m= n / ps[i];
d=extEuclid(ps[i], m, x, y);
ret=(ret+y*m*as[i])%n;
}
return (n+ret%n)%n;
}
ll gs[55],Ps[55],ss[55],piR[55],aiR[55],pcntR;
ll ans[100005],nans,tans[100005],tnans;// error1
void merge(ll sol,ll prevPi, ll curPi){
as[0] = sol, ps[0] = curPi;
for (int j = 1; j <= nans; ++j) {
as[1] = ans[j], ps[1] = prevPi;
tans[++tnans] = crt(as,ps,2);
}
}
void kResidue(ll k, ll y, ll P){
nans=0;
ll pi2_ai=1;
if(!(P & 1ll)){
pi2_ai=2;
P>>=1;
if(y&1) ans[++nans]=1;
else ans[++nans]=0;
while(!(P&1ll)){
pi2_ai*=2;
tnans=0;
for(int i=1;i<=nans;++i){
if(qpow(ans[i], k, pi2_ai) == y % pi2_ai) tans[++tnans]=ans[i];
if(qpow(ans[i]+(pi2_ai >> 1), k, pi2_ai) == y % pi2_ai) tans[++tnans]= ans[i] + (pi2_ai >> 1);
}
for(int i=1;i<=tnans;++i)
ans[i]=tans[i];
if(tnans==0) {
ans[nans = 1] = -1;return; //error5
}
nans=tnans;
P>>=1;
}
}
factor(P, piR, aiR, pcntR);
ll privPi=pi2_ai;//error
for(int i=1; i <= pcntR; ++i){
tnans = 0;
ll pi_ai=pow(piR[i],aiR[i]);
ll divisor=gcd(y%pi_ai,pi_ai);
if(divisor==pi_ai){
ll t = upperInt(aiR[i],k);
ll base=pow(piR[i],t);
for(ll sol=base;sol<=pi_ai;sol+=base){// error6 base==pi_ai, error7 +=base
if(i==1 && pi2_ai == 1){
tans[++tnans]=sol;
}else {
merge(sol,privPi,pi_ai);
}
}
}else{
ll phi=pi_ai/piR[i]*(piR[i]-1);
ll upperphi=phi;
ll yy=y%pi_ai;
ll tpi_ai=pi_ai;
ll pek=1;
if(divisor!=1){
int tmp=0;
while((divisor/pek)%piR[i]==0){
pek*=piR[i];
tmp++;
}
if(tmp/k*k!=tmp){
ans[nans=1]=-1;return ;
}
tpi_ai/=pek;
phi=tpi_ai/piR[i]*(piR[i]-1);
yy/=pek;
pek=qpow(piR[i],tmp/k,pi_ai);
}
ll g=getG(tpi_ai,phi);
ll c=bsgs(qpow(g,k,tpi_ai),yy,tpi_ai);
if(c==-1){
ans[nans=1]=-1;return ;
}
ll delta= phi / gcd(k, phi);
for(ll cur=c%delta,sol; cur < upperphi; cur+=delta){
sol=(qpow(g, cur, pi_ai)*pek)%pi_ai;
if(i==1 && pi2_ai == 1){
tans[++tnans]=sol;
}else{
merge(sol,privPi,pi_ai);
}
}
}
sort(tans+1,tans+1+tnans);
tnans=unique(tans+1,tans+1+tnans)-tans-1;
for(int j=1;j<=tnans;++j){
ans[j]=tans[j];
}
nans=tnans;
privPi*=pi_ai;//error4
}
sort(ans+1,ans+1+nans);
nans=unique(ans+1,ans+1+nans)-ans-1;
return ;
}
int main() {
// int v7[4]; // [rsp+10h] [rbp-50h]
// int v8[4]; // [rsp+20h] [rbp-40h]
// int v9[4]; // [rsp+30h] [rbp-30h]
// v7[0] = -855924449;
// v7[1] = 2093860801;
// v7[2] = 767425486;
// v7[3] = 480046732;
// v8[0] = 873063593;
// v8[1] = 947071579;
// v8[2] = 1040105093;
// v8[3] = 780771613;
// v9[0] = -785653175;
// v9[1] = -1371619377;
// v9[2] = -1173543495;
// v9[3] = -1099902543;
// for(int _=0;_<4;++_){
// kResidue((unsigned int)v8[_],(unsigned int)v7[_],(unsigned int)v9[_]);
// if(ans[1]==-1){
// puts("No Solution");
// continue;
// }
// for(int i=1;i<=nans;++i){
// printf("%lld ",ans[i]);
// }
// puts("");
// }
ll k,y,M;
int cs=0;
scanf("%d",&cs);
for(int _=0;_<cs;++_){
scanf("%lld %lld %lld",&k,&M,&y);
kResidue(k,y,M);
if(ans[1]==-1){
puts("No Solution");
continue;
}
for(int i=1;i<=nans;++i){
printf("%lld ",ans[i]);
}
puts("");
}
return 0;
}
/*
1
90 100 0
1
3 36 17
Input
1
3 531441 330750
Output
264 19947 39630 59313 78996 98679 118362 138045 157728 177411 197094 216777 236460 256143 275826 295509 315192 334875 354558 374241 393924 413607 433290 452973 472656 492339 512022
*/