您的位置:首页 > 其它

高次同余方程模板BabyStep-GiantStep

2015-12-07 00:23 417 查看
/*************************************
---高次同余方程模板BabyStep-GiantStep---

输入:对于方程A^x=B(mod C),调用BabyStep(A,B,C),(0<=A,B,C<=10^9)

输出:无解放回-1,有解放回最小非负整数x

复杂度:O(C^0.5),只与C有关,与A,B的大小无关
************************************/

typedef long long ll;
#define HASH_N 100007

struct hashnode
{
int next;
ll key;
int j;
}HashLink[ HASH_N ];

int hashpre[ HASH_N ],hashcnt;

void hash_insert(ll x,ll key,int j)
{
for(int p=hashpre[x];p!=-1;p=HashLink[p].next)
{
if(HashLink[p].key==key) return ;
}
HashLink[ hashcnt ].key = key;
HashLink[ hashcnt ].j = j;
HashLink[ hashcnt ].next = hashpre[x];
hashpre[x] = hashcnt++;
}

int hash_get(ll key)
{
ll x=key%HASH_N;
for(int p=hashpre[x];p!=-1;p=HashLink[p].next)
{
if( HashLink[p].key == key ) return HashLink[p].j;
}
return -1;
}

ll gcd(ll a,ll b)
{
if(b==0) return a;
return gcd(b,a%b);
}

//ax + by = gcd(a,b)
//传入固定值a,b.放回 d=gcd(a,b), x , y
void extendgcd(long long a,long long b,long long &d,long long &x,long long &y)
{
if(b==0){d=a;x=1;y=0;return;}
extendgcd(b,a%b,d,y,x);
y-=x*(a/b);
}

//Ax=1(mod M)
//返回x的范围是[0,M-1]
ll GetNi(ll A,ll M)
{
ll rex=0,rey=0;
ll td=0;
extendgcd(A,M,td,rex,rey);
return (rex%M+M)%M;
}

//a^b%mod 快速幂
long long Quk_Mul(long long a,long long b,long long mod)
{
long long qsum=1;
while(b)
{
if(b&1) qsum=(qsum*a)%mod;
b>>=1;
a=(a*a)%mod;
}
return qsum;
}

//测试x较小的情况,必须!
ll firsttest(ll A,ll B,ll C)
{
ll tmp=1;
if(B==1) return 0;
for(int i=1;i<100;i++)
{
tmp = (tmp*A)%C;
if(tmp==B) return i;
}
return -1;
}

//欧拉函数 返回[1,x-1]中与x互素的数的个数,复杂度n^0.5
ll euler(ll x)
{
ll i, res=x;
for (i = 2; i < (ll)sqrt(x * 1.0) + 1; i++)
if(x%i==0)
{
res = res / i * (i - 1);
while (x % i == 0) x /= i;
}
if (x > 1) res = res / x * (x - 1);
return res;
}

ll BabyStep(ll A,ll B,ll C)
{
if(0==A || 0==C) return -1;
if(C==1) return 0;
B = B%C;
ll ans = firsttest(A,B,C);//为了防止x比较小的时候
if(ans != -1) return ans;
ll D=1;
int c=0;
ll d;
while( (d=gcd(A,C)) != 1 )
{
if( B%d !=0 ) return -1;//无解的情况
B /= d;
C /= d;
D = D*A/d%C;
c++;
}

//得到了 D*A^(x-c)=B (mod C) ,gcd(A,C)=1 , gcd(D,C)=1
ll D_1=GetNi(D,C);//求D的逆元
B = B*D_1%C;
//求A^x=B (mod C),然后返回x+c
ll m = ceil( sqrt(C+0.0) );

memset(hashpre,-1,sizeof(hashpre));
hashcnt=0;
ll hashnum=1;
hash_insert(1, 1, 0);
for(int i=1;i<m;i++)
{
hashnum = (hashnum*A)%C;
hash_insert(hashnum%HASH_N, hashnum ,i);
}

ll ol=euler(C);
ll mA=Quk_Mul(A, m, C);
ll ta=1;

ll tmpni = Quk_Mul(mA, ol-1, C);

for(int i=0;i<m;i++,ta=ta*tmpni%C)
{
//解线性同余方程 tx=B*ta (%C) ,ta直接用费马小定理求逆元
ll tx = ta;
tx = (tx*B)%C;
int j=hash_get(tx);
if(j!=-1)//找到解了
{
return i*m+j+c;
}
}

return -1;
}


其实还是很很正常的解法,这个甚至可以当成一道难一点的数论题目做。

关于详细的解释。

[转自http://blog.csdn.net/auto_ac/article/details/11842695]

题意:求满足a^x=b(mod n)的最小的整数x。

分析:很多地方写到n是素数的时候可以用Baby step,Giant step, 其实研究过Baby step,Giant step算法以后,你会发现 它能解决 “n与a互质”的情况,而并不是单纯的n是素数的情况。如果a与n不是互质的,那么我们需要处理一下原方程,让a与n互质,然后再用Baby step,Giant step解出x即可。

Baby step,Giant step算法思想:对于a与n互质,那么则有a^phi(n)=1(mod n), 对于n是素数phi(n) == n-1, 否则phi(n) < n-1, 所以x的取值只要在0----n-2之中取就可以了。

当n很小时,可以直接枚举,但当n很大时,肯定会超时,Baby step,Giant step就是用了一种O(sqrt(n)*log(n))的方法枚举了所有的0-----n-2。令m = sqrt(n);

我们可以预处理出a^0,a^1,.........a^m,都放入哈希表中, 然后 (a^m)^i+v(哈希表里的其中一个值)就一定是解,每次枚举i(0-----m-1),计算出v,判断v是否出现在哈希表中,如果有就是解。 对于m为什么取sqrt(n)是为了复杂度的平衡,这一点是跟分块算法很相似的。

对于a与n不互质的情况分析:令 t = gcd(a,n),那么a与n都约去t,当然b也要约去t(不能约去就无解),约去一个t以后方程就变为 aa*a^(x-1) = bb(mod nn), (其中 aa = a/t bb = b/t nn = n/t) , 这里nn还可能与a不互质,那么我们一直拿出一个新的a对(a, bb, nn)约去t,直到a与nnn....(nnn...表示约去若干次t以后的n)互质。以下用(用三个字母表示约去若干次后,如bbb) 则结果为aa^c*a^(x-c) = bbb(mod nnn), 我们让等式左右分别乘以aa^c关于nnn的逆元 变为a^(x-c) = w (mod nnn) , w =bbb *(aa^c)^(-1)。 a^x = w (mod n)可以用bbb *(aa^c)^(-1)Baby step,Giant step直接求出,如果有解那把未知数+c。

具体看代码中的cal函数。

注意:在以上过程中x有可能<c,所以我们必须每约去一个t就要特判一下当前情况aa 与 bb就说明当前c是解。

哈希表实现看题目时间要求,map太慢,自己手写hash是很快的。
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: