高次不定方程BSGS算法
学习数学真是一件赛艇的事.
BSGS名字听起来非常有意思,力拔山兮气盖世,北上广深,小步大步...算法其实更有意思,它是用来求解一个方程的
\(A^x≡B mod P\)
是不是特别眼熟,有几个式子长的特别像,先观察一下:
一:快速幂: 求\(A^B mod P\)的值
二:乘法逆元 \(A*x ≡ 1 (mod P)\)
或者 \(A*x ≡ B (mod P)\)
三:欧拉定理 \(A^{φ(P)}≡ 1 (mod P)\) (A,P互质)
四:费马小定理 \(A^{P-1} ≡ 1 (mod P)\) (P是质数)
先说下这四者关系:快速幂可以快速求后三个,费马小定理是欧拉定理的特殊情况,逆元可以通过费马小定理和快速幂解
可如果有这样的一个式子:
\(A^x ≡ B (mod P)\) 我们先假设A,P互质
好像和这四个式子都很像,所以呢?所以呢?
当年我们证明费马小定理的时候发现这个x的范围是在\([0,p-1]\)之间
那么我们就可以枚举x从0到p-1,复杂度为O(P);
应该能拿上30分
接下来就是一波骚操作:我们令 \(m=\lceil\sqrt{P}\rceil\) ,然后就可以设 \(x = i*m+j\) ,其中\(i=\lfloor \frac{x}{m} \rfloor ,j=x%m\),把x代入原来的式子可以得到$A ^{i*m+j} ≡ B ( mod P ) $两边乘上一个 \(A^{-i*m}\),就可以得到 \(A^j≡B*A^{-i*m}( mod P )\)
所以呢?
所以就可以求了啊
我们只要枚举左边的\(j\),把左边的答案和\(j\)存起来 \(left\)_\(ans [ j ] = A ^ j % P\) (可是存不下怎么办,哈希蛤一下就存下了)然后再枚举右边的 i,计算右边的值,看看我们右边的值是否在数组里出现过,如果出现过那么我们通过i和j找到的 i*m+j 就是一个答案了
然后就会发现复杂度被我们开了一个方
冷静分析:
这个算法先枚举j需要\(\sqrt{P}\)的时间,再枚举i需要\(\sqrt{P}\)的时间,不过枚举i是要算下逆元需要\(log_2P\)的时间,看起来复杂度=O(\(\sqrt{P}\)+\(\sqrt{P}\)log2(P))=O(\(\sqrt{P}\)\(log_2P\)), 不过我们再看看右边的式子:
$ B A^{-im} mod P = B *(A{-m})imod P $
然后我们就得到右边的递推式,只要先求出 \(A^{-m} modP\) 就可以O(1)计算右边的式子了,其中\(A^{-m}≡A^{P-1-m}modP\),因为费马小定理...所以复杂度被我们降到了O(\(\sqrt{P}\)+\(log_2P\)+\(\sqrt{P}\))=O(\(\sqrt{P}\))
灼热分析:
算法的思想其实就是分块,把x分成\(\sqrt{P}\)\(\sqrt{P}\)的块,会到设x的式子,x=im+j,我们先Baby_Step枚举小的j,再Giant_Step枚举大的i,名字听起来很形象哈哈哈哈哈哈.因为先枚举小的,后枚举大的,所以当出现i满足条件时,可以保证此时答案是最小的正整数解,这时直接return i*m+j
科学分析:
哈希好用呐~之前懒得用哈希,总觉得用STL的map能省很多事,然后就很尴尬的调了两天...一直TLE,最后绝望的手写了哈希表,然后居然就p+的A掉了,千万别用map,千万别用map,千万别用map,STL里面的玄学操作看起来很好用,我们最好还是乖乖学一学正常操作,老老实实手写哈希....
然后看看代码
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>
using namespace std;
#define ll long long
const int mod=1048573;
int hcnt=0,head[mod+10];
struct Haha{
int val,id,next;
}hash[mod+10];
void insert(int x,int pos){
int k=x%mod;
hash[++hcnt].val=x;
hash[hcnt].id=pos;
hash[hcnt].next=head[k];
head[k]=hcnt;
}
int find(int x){
int k=x%mod;
for(int i=head[k];i;i=hash[i].next){
if(hash[i].val==x) return hash[i].id;
}
return -1;
}
ll ksm(int a,int b,int p){
int x=a;
ll ret=1;
if(b<0) return -1;
while(b){
if(b&1) ret=1ll*(ret*x)%p;
b>>=1;
x=1ll*x*x%p;
}
return ret;
}
int BSGS(int a,int b,int p){
int m=(int)(sqrt(p)+0.999999);
if(b==1) return 0;
if(a==b) return 1;
if(!b){
if(!a) return 1;
return -1;
}
ll x=1;
for(int i=1;i<=m;++i){
x=x*a%p;
insert(x,i);
}
ll inv=1;
int inv2=ksm(a,p-m-1,p)%p;
for(int i=0;i<m;++i){
int k=i*m;
if(inv==-1) return -1;
int ans=inv*b%p;
int jgy=find(ans);
if(~jgy){
return k+jgy;
}
inv=1ll*inv*inv2%p;
}
return -1;
}
void init(){
for(int i=0;i<mod;++i){
hash[i].val=-1;
hash[i].next=0;
hash[i].id=0;
}
memset(head,0,sizeof(head));
hcnt=0;
}
int main(){
int a,b,p;
while(~scanf("%d%d%d",&p,&a,&b)){
init();
int dove=BSGS(a,b,p);
if(~dove) printf("%d\n",dove);
else printf("no solution\n");
}
return 0;
}
