参考:http://blog.miskcoo.com/2015/05/discrete-logarithm-problem
大步小步算法。
求模方程
ax≡b(modp)
的最小解。p为素数。
分块思想。
popoqqq说是meet in the middle。有时间学学那个东西吧。
我们要找的
x=A⌈p√⌉+B
模方程变成:
aA⌈p√⌉+B≡b(modp)
不难看出,如果我们预处理出前
B
个
以后可以每
⌈p√⌉
个值里,模方程左右同时乘上
a−⌈p√⌉
。
t次后,模方程变为
aB≡b∗at∗−⌈p√⌉(modp)
。
只需要找我们预处理出的前B个值里是否有
b∗at∗−⌈p√⌉
即可。用哈希表做。
(因为p是素数所以
a−⌈p√⌉
一定存在)
若第一次预处理出
m
个值,那么时间复杂度为
bzoj3239
/**************************************************************
Problem: 3239
User: flyfish0851
Language: C++
Result: Accepted
Time:212 ms
Memory:2404 kb
****************************************************************/
#include<cstdio>
#include<cmath>
#include<map>
#define ll long long
using namespace std;
map<ll,ll>hash;
ll pow(ll a,ll b,ll p){
ll ret=1;
while(b){
if(b&1) ret=ret*a%p;
a=a*a%p;b>>=1;
}
return ret;
}
ll inv(ll a,ll p){return pow(a,p-2,p);}
ll BSGS(ll a,ll b,ll p){
hash.clear();
ll m=(int)sqrt(p)+1;
ll e=1;hash[e]=0;
ll v=inv(pow(a,m,p),p);
for(int i=1;i<m;i++){
e=e*a%p;
if(!hash.count(e)) hash[e]=i;
}
for(int i=0;i<=m;i++){
if(hash.count(b)) return (ll)i*m+hash[b];
b=b*v%p;
}
return -1;
}
int main(){
//freopen("in.txt","r",stdin);
ll a,b,p;
while(~scanf("%lld%lld%lld",&p,&a,&b)){
ll t=BSGS(a,b,p);
if(t==-1) puts("no solution");
else printf("%lld\n",t);
}
return 0;
}
扩展大步小步。
求模方程
ax≡b(modm)
的最小解。m任意。
这次少了m为素数的限制。
我们发现瓶颈在求
a−⌈p√⌉
,它不一定存在。我们的做法是把模方程同时除以
g=gcd(a,p)
m′=mg,b′=bg(ag)−1
新的模方程为
ax′≡b′(modm′)
最终的
x=x′+1
(此过程中a不变。)
不断重复这个过程最后一定会得到一个可以解的方程,套用刚刚的大步小步法解出后即可。要注意的是在这个过程中如果某一步发现
b′=1
,那么就可以直接退出,因为这时候已经得到了解。
代码有几个需要注意的细节:
1、
a=0
特判
2、
p=1
3、
b=b∗inv(a/g,p)mod p
虽然我说的细节对于很多人都不是细节吧= =
bzoj2480&bzoj1467
#include<cstdio>
#include<cmath>
#include<map>
#define ll long long
using namespace std;
map<ll,ll>hash;
ll gcd(ll a,ll b){return b==0?a:gcd(b,a%b);}
void exgcd(ll a,ll b,ll &x,ll &y){
if(!b) {x=1;y=0;return;}
exgcd(b,a%b,y,x);y-=x*(a/b);
}
ll inv(ll a,ll p){if(!a) return 1;ll x,y;exgcd(a,p,x,y);return (x%p+p)%p;}
ll pow(ll a,ll b,ll p){
ll ret=1;
while(b){
if(b&1) ret=ret*a%p;
a=a*a%p;b>>=1;
}
return ret;
}
ll BSGS(ll a,ll b,ll p){
if(a==0) return b==0?1:-1;
b%=p;a%=p;
hash.clear();
ll cnt=0;
for(int g=gcd(a,p);g!=1;g=gcd(a,p)){
if(b%g) return -1;
b/=g;p/=g;b=b*inv(a/g,p)%p;
cnt++;
if(b==1) return cnt;
}
if(p==1) return cnt;
int m=(int)sqrt(p)+1;
ll e=1;hash[e]=0;
ll v=inv(pow(a,m,p),p);
for(int i=1;i<m;i++){
e=e*a%p;
if(!hash.count(e))
hash[e]=i;
}
for(int i=0;i<=m;i++){
if(hash.count(b)) return i*m+hash[b]+cnt;
b=b*v%p;
}
return -1;
}
int main(){
ll a,b,p;
while(scanf("%lld%lld%lld",&a,&p,&b)&&p){
ll t=BSGS(a,b,p);
if(t==-1) puts("No Solution");
else printf("%lld\n",t);
}
return 0;
}