离散对数之BSGS(扩展BSGS)与Pohlig-Hellman(知识点整理+板子总结)

最新推荐文章于 2022-07-20 13:53:39 发布

Code92007

最新推荐文章于 2022-07-20 13:53:39 发布

阅读量1k

点赞数 2

分类专栏：知识点总结文章标签：知识点总结离散对数 BSGS Pohlig-Hellman CRT

本文链接：https://blog.csdn.net/Code92007/article/details/98661838

版权

知识点总结专栏收录该内容

89 篇文章 6 订阅

订阅专栏

心得

离散对数这里还是不怎么熟，

很早就学过，然而不怎么会套板子，

现在总算是把板子整理下来了

思路来源

https://blog.csdn.net/qq_34921856/article/details/79794335

https://blog.csdn.net/litble/article/details/73252928

知识点整理

1.BSGS

求 $a^{x}\equiv b(mod\ p)$ 的最小x值，p仅为素数

令 x=i*m+j,j<m ，则需解 $a^{i*m}\equiv b*a^{-j}(mod\ p)$ ，

j从0到m-1预处理等式右边的值v，令map[v]=j存下标，

再左端枚举i从0到m-1，去查map表即可，找到最小的i满足j在表中，i*m+j即为答案

复杂度最低时，令 $m=\sqrt p$ ， $O(\sqrt p)$ 解决一个x

/*==================================================*\
| Baby-Step-Giant-Step 大步小步算法
| 求 a^x === b (mod p) 中的 最小 x值 -- 此处p仅为素数
| 实际运用中用自己的hash表代替map可防TLE
\*==================================================*/
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
ll a,b,p;
ll extgcd(ll a,ll b,ll &x,ll &y)
{
	ll d=a;
	if(b)d=extgcd(b,a%b,y,x),y-=(a/b)*x;
	else x=1,y=0;
	return d;
}
ll BSGS(ll a, ll b, ll p) 
{ 
    a %= p; b %= p;
    map<ll, ll> h;
    ll m = ceil(sqrt(p)), x, y, d, t = 1, v = 1;
    for(ll i = 0; i < m; ++i)
	{
        if(h.count(t)) h[t] = min(h[t], i);
        else h[t] = i;
        t = (t*a) % p;
    }
    for(ll i = 0; i < m; ++i) 
	{
        d = extgcd(v, p, x, y);//vx==1(mod p)即v在模p意义下逆元 
        x = (x* b/d % p + p) % (p);
        if(h.count(x)) return i*m + h[x];
        v = (v*t) % p;
    }
    return -1;
}
int main()
{
	scanf("%lld%lld%lld",&a,&b,&p);//a^x==b(mod p)
	printf("%lld\n",BSGS(a,b,p));
	return 0;
}

2.扩展BSGS

求 $a^{x}\equiv b(mod\ p)$ 的最小x值，p不为素数

$a^{x-1}*a\equiv b(mod\ p)$ ，

考虑扩展欧几里得，在求 $ax\equiv b(mod\ p)$ 时，如果a、p不互质则a、b、p同除gcd(a,p)

以上方法也同理， $a^{x-1}*\frac{a}{d}{}\equiv \frac{b}{d}(mod\ \frac{p}{d})$ ，

若 $\frac{p}{d}$ 与a互质，直接解x-1，所得加1即答案

否则再提一个a出来， $a^{x-2}*\frac{a}{d_{1}}*\frac{a}{d_{2}}\equiv \frac{b}{d_{1}*d_{2}}(mod\ \frac{p}{d_{1}*d_{2}})$ ，

重复该过程，直至a与 $\frac{p}{\prod_{i=1}^{k}d_{i}}$ 互质，

去解 $a^{x-k}*\frac{a^{k}}{\prod_{i=1}^{k}d_{i}} \equiv \frac{b}{\prod_{i=1}^{k}d_{i}}(mod\ \frac{p}{\prod_{i=1}^{k}d_{i}})$

注意解出来的是x-k，再加k即答案，

这里整理的板子，是手写哈希的，开大于预处理部分的空间的最小素数值

#include<iostream>
#include<cstdio>
#include<climits>
#include<algorithm>
#include<cstring>
#include<cmath>
using namespace std;
#define ll long long
#define mod 99991
int head[mod],tot;
ll Next[mod<<1],id[mod<<1],has[mod<<1];
ll gcd(ll x,ll y)
{
	ll r=x%y;
	while(r){x=y;y=r;r=x%y;}
	return y;
}
void add(ll x,ll num)
{
	ll k=num%mod;
	id[++tot]=x;has[tot]=num;
	Next[tot]=head[k];head[k]=tot;
}
ll find(ll num)
{
	ll k=num%mod;
	for(ll i=head[k];i!=-1;i=Next[i])
	if(has[i]==num)return id[i];
	return -1;
}
ll exbsgs(ll a,ll b,ll p)
{
	a%=p;b%=p;
	if(b==1)return 0;
	ll tmp=1,d=1,cnt=0;
	while((tmp=gcd(a,p))!=1)
	{
		if(b%tmp)return -1;
		cnt++;b/=tmp;p/=tmp;d=d*(a/tmp)%p;
		if(b==d)return cnt;//注意这里
	}
	ll m=ceil(sqrt(p)),q=1,j;
	//一般m开为根号p 视具体询问而改变 
	//看预处理和询问的次数 如只询问一次则应根号p 
	for(ll i=0;i<m;i++)
	{add(i,(q*b)%p);q=(q*a)%p;}
	for(ll i=1;i<=m+1;i++)
	{
		d=(d*q)%p;j=find(d);
		if(j!=-1)return i*m-j+cnt;
	}
	return -1;
}
void init()
{
	memset(head,-1,sizeof head);
	tot=0;
}
int main()
{
	ll a,b,p,ans;
	while(~scanf("%lld%lld%lld",&a,&p,&b))
	{
		if(!a&&!b&&!p)break;
		init();
		ans=exbsgs(a,b,p);
		if(ans==-1)puts("No Solution");
		else printf("%lld\n",ans);
	}
	return 0;
}

3.Pohlig-Hellman

http://www-math.ucdenver.edu/~wcherowi/courses/m5410/phexam.html

啃了一天，贴一下网页局部，描述一下算法流程，讲点心得

求 $a^{x}\equiv b(mod\ p)$ 的最小x值，p仅为素数但p<=1e18，此时 $O(\sqrt p)$ 解决不了x，

但p的特点在于，p-1可以表示为小质数的乘积，如例 $p-1=2^{2}*3^{4}*5^{2}$ ，

如果第一步能先求出x模每个素因子的幂 $p_{i}^{k_{i}}$ 的值，届时第二步CRT合并即可

考虑第一步，

要求最小循环节，初始循环节cycle=p-1

应当检查每个素因子是否满足，对于当前循环节cycle而言，

$a^{\frac{cycle}{p_{i}}}\equiv 1(mod\ p)$ 是否成立，如果成立，

说明cycle可以除掉 $p_{i}$ ，使循环节更小，

那就除掉 $p_{i}$ ，直至该素因子 $p_{i}$ 不满足 $a^{\frac{cycle}{p_{i}}}\equiv 1(mod\ p)$ ，

再去检查下一素因子，直至检查完所有素因子

经验证，这里的p的最小循环节为p-1，

先预处理 $a^{1*\frac {p-1}{p_{i}}},...a^{(p_{i}-1)*\frac {p-1}{p_{i}}}$ ，对于x2而言，只有 $a^{1*\frac{p-1}{2}}\equiv 8100(mod\ 8101)$ ，

然后去求 $(a^{x})^{\frac{p-1}{2}}\equiv7531^{\frac{p-1}{2}}\equiv8100(mod\ p)$ ，

$x_{2}$ 是x除以4的余数，可以表示为 $x_{2}=c_{0}+c_{1}*2$

即 $x=c_{0}+c_{1}*2+c_{2}*4+...$ ，显然后面包含*2的项是循环节的倍数，被循环节消掉了，

有 $(a^{x})^{\frac{p-1}{2}}\equiv(a^{c_{0}+c_{1}*2+c_{2}*4+...})^{\frac{p-1}{2}}\equiv a^{c_0*\frac{p-1}{2}}$ ，查表知 $c_{0}=1$

对原余数7531乘以 $a^{c_{0}}$ 的逆元，新的 $x=c_{1}*2+c_{2}*4+...$ ，重复该过程，

$(a^{x})^{\frac{p-1}{4}}\equiv(a^{c_{1}*2+c_{2}*4+...})^{\frac{p-1}{4}}\equiv a^{c_1*2*\frac{p-1}{4}}\equiv a^{c_1*\frac{p-1}{2}}\equiv 1$ ，查表知 $c_{1}=0$

由于 $x_{2}=c_{0}+c_{1}*2$ ，代入解得 $x_{2}=1$ ，该过程终止

对于每个素因子 $p_{i}^{k_{i}}$ ，都类似地解出其余数 $x_{i}$ ，得到同余方程组 $\begin{bmatrix} x \equiv x_{2}(mod\ 2^{2})\\ x \equiv x_{3}(mod\ 3^{4}) \\ x \equiv x_{5}(mod\ 5^{2}) \end{bmatrix}$ ，

由于素因子幂次之间显然互质，裸的CRT，合并构造一下即得x，

附：hdu6632 discrete logarithm problem

求最小的x满足 $a^{x}\equiv b(mod\ p)$ ，p-1的素因子只由2、3组成

$2<=a,b<=p-1,65537<=p<=10^{18}$

~~据说是CTF界人尽皆知sb题，然而我一个打ACM的怎么可能会知道~~

~~（抓周树人和我鲁迅有什么关系）~~

仿着网页的过程敲，写了一波又丑常数又大的代码

#include<bits/stdc++.h>
using namespace std;
const int maxn=1e4+10;
const int N=60;
const int MAXP=20;//最大素因子的大小 
typedef long long ll;
bool ok[maxn],yes;
int prime[maxn],cnt;
int t;
int v[N],num[N],tot;//素数因子v[i] 出现了num[i]次 
ll X[N],P[N],sz;//返回(x',pi^ki) x'=x(mod pi^ki) 用于CRT 
ll a,b,p;
void sieve()
{
	for(ll i=2;i<maxn;++i)
	{
		if(!ok[i])prime[cnt++]=i;
		for(int j=0;j<cnt;++j)
		{
			if(i*prime[j]>=maxn)break;
			ok[i*prime[j]]=1;
			if(i%prime[j]==0)break; 
		}
	}
}
int fac(ll x)
{
	int tot=0;
	for(int i=0;i<cnt;++i)
	{
		if(1ll*prime[i]*prime[i]>x)break;
		if(x%prime[i]==0)
		{
			v[++tot]=prime[i];
			num[tot]=0;
			while(x%prime[i]==0)
			{
				num[tot]++;
				x/=prime[i];
			}
		}
	}
	if(x>1)
	{
		v[++tot]=x;
		num[tot]=1;
	}
	return tot;
}
ll extgcd(ll a,ll b,ll &x,ll &y)
{
	ll d=a;
	if(b)d=extgcd(b,a%b,y,x),y-=(a/b)*x;
	else x=1,y=0;
	return d;
}
ll mul(ll u,ll v,ll p)
{
	return (u*v-(ll)((long double)u*v/p)*p+p)%p;
}
ll modpow(ll x,ll n,ll mod)
{
	ll res=1;
	for(;n;n/=2,x=mul(x,x,mod)%mod)
	if(n&1)res=mul(res,x,mod)%mod;
	return res;
}
ll Pholig_Hellman(ll a,ll B,ll p)//返回(x',pi^ki) x'=x(mod pi^ki) 
{
	ll r=p-1;
	//找a的阶(最小)r p-1不一定是最小 类似判原根的过程 
	//最后r应该化为最小的阶
	tot=fac(r);
	for(int i=1;i<=tot;++i)
	while(r%v[i]==0&&modpow(a,r/v[i],p)==1)r/=v[i];
	tot=fac(r);//重新求素因子 
	ll x,y,inv,tmp[MAXP],sz=0;
	extgcd(a,p,x,y);
	inv=(x%p+p)%p;
	for(int i=1;i<=tot;++i)
	{
		int c=v[i],d=num[i];
		ll b=B,f=0,g=1;
		ll l=1,now=r,h=modpow(a,now/c,p);
		//printf("h:%lld\n",h);
		for(int k=0;k<c;++k,l=mul(l,h,p))
		{
			tmp[k]=l; 
			//printf("tmp[%d]:%lld\n",k,tmp[k]);
		}
		for(int j=1;j<=d;++j)
		{
			now/=c;//(p-1)/(c^j)
			ll e=modpow(b,now,p);
			ll co=-1,dif;//系数co 
			for(int k=0;k<c;++k)
			{
				if(e==tmp[k])
				{
					co=k; 
					break;
				}
			}
			if(co==-1)return -1;//无解
			//printf("co:%lld\n",co);
			co*=g;//实际系数 
			f+=co;
			dif=modpow(inv,co,p);
			b=mul(b,dif,p); 
			g*=c;
		}
		//printf("f:%lld g:%lld\n",f,g);
		X[++sz]=f;
		P[sz]=g;
	}
	return sz;
}
ll CRT(ll r[],ll m[],int n)//x'==x(mod p)方程组 共n项 要求pi互素 
{
	ll ans=0,M=1,x,y;
	for(int i=1;i<=n;++i)
	M*=m[i];
	for(int i=1;i<=n;++i)
	{
		ll mi=M/m[i];
		extgcd(mi,m[i],x,y);
		x=(x%m[i]+m[i])%m[i];
		ans=(ans+mul(mul(r[i],mi,M),x,M))%M;
	} 
	return (ans+M)%M;
}
int main()
{
	sieve();
	scanf("%d",&t);
	while(t--)
	{
		scanf("%lld%lld%lld",&p,&a,&b); 
		sz=Pholig_Hellman(a,b,p);
		if(sz!=-1)printf("%lld\n",CRT(X,P,sz));
		else puts("-1");
	}
	return 0;
}
//p a b
/*
105
65537 2 4
8101 6 7531
*/

Code92007

关注

2
点赞
踩
9

收藏

觉得还不错? 一键收藏
打赏
0
评论
离散对数之BSGS(扩展BSGS)与Pohlig-Hellman(知识点整理+板子总结)

心得离散对数这里还是不怎么熟，很早就学过，然而不怎么会套板子，现在总算是把板子整理下来了思路来源https://blog.csdn.net/qq_34921856/article/details/79794335https://blog.csdn.net/litble/article/details/73252928知识点整理1.BSGS求的最小x值，p仅为素数...
复制链接

扫一扫