咳咳,玩了一回标题党哈~
这神奇的题目就是——
codevs1565/洛谷P2485/bzoj2242【SDOI2011】计算器
#快速幂
这个…真的需要讲吗?已经在学扩展欧几里德和bsgs的oier应该是早就会快速幂的吧…
其原理如下。
若y为偶数,则有
x
y
=
(
x
y
/
2
)
2
x^y=(x^{y/2})^2
xy=(xy/2)2
若y为奇数,则有
x
y
=
(
x
y
/
2
)
2
∗
x
x^y=(x^{y/2})^2*x
xy=(xy/2)2∗x (y/2是向下取整)
然后就得到了写法:
ll ksm(ll y,ll z,ll p){//求(y^z)%p
ll re=1;y%=p;
while(z){
if(z&1)re=(y*re)%p;
z>>=1;y=(y*y)%p;
}
return re%p;
}
#扩展欧几里德
这又是一种神奇的算法
我们知道,
a
x
≡
b
(
ax≡b(
ax≡b(mod
c
)
c)
c)可以化为方程
a
x
+
c
y
=
b
ax+cy=b
ax+cy=b,这个方程有整数解的前提是
g
c
d
(
a
,
c
)
∣
b
gcd(a,c)|b
gcd(a,c)∣b,然后我们把这个神奇的方程解出
x
x
x最小非负整数解即可。
怎么解?用扩展欧几里德啊!
扩展欧几里德可以用来算
a
x
+
c
y
=
g
c
d
(
a
,
c
)
ax+cy=gcd(a,c)
ax+cy=gcd(a,c),因为
g
c
d
(
a
,
b
)
=
g
c
d
(
b
,
a
gcd(a,b)=gcd(b,a
gcd(a,b)=gcd(b,a%
b
)
b)
b)所以
p
∗
a
+
q
∗
b
=
g
c
d
(
a
,
b
)
=
g
c
d
(
b
,
a
p*a+q*b=gcd(a,b)=gcd(b,a
p∗a+q∗b=gcd(a,b)=gcd(b,a%
b
)
=
p
∗
b
+
q
∗
(
a
b)=p*b+q*(a
b)=p∗b+q∗(a%
b
)
=
p
∗
b
+
q
(
a
−
a
/
b
∗
b
)
=
q
∗
a
+
(
p
−
a
/
b
∗
q
)
∗
b
b)=p*b+q(a-a/b*b)=q*a+(p-a/b*q)*b
b)=p∗b+q(a−a/b∗b)=q∗a+(p−a/b∗q)∗b这样我们就可以一边求gcd一遍玩扩展欧几里德啦!
ll exgcd(ll a,ll b,ll &r1,ll &r2){//a*r1+b*r2=gcd(a,b)
if(!b){r1=1;r2=0;return a;}
ll re=exgcd(b,a%b,r1,r2),tmp;
tmp=r1;r1=r2;r2=tmp-(a/b)*r2;
return re;
}
玩出来
x
x
x了之后,再让
x
=
x
∗
(
b
/
g
c
d
(
a
,
c
)
)
x=x*(b/gcd(a,c))
x=x∗(b/gcd(a,c))就可以啦!
#bsgs
bsgs算法,又称为大步小步算法,拔山盖世算法,北上广深算法,百事公司算法等等…
这个算法是这样的:
求使得
y
x
≡
z
(
y^x≡z(
yx≡z(mod
p
)
p)
p)的最小
x
x
x,(
p
p
p是质数),由费马小定理可知,
y
p
−
1
≡
1
(
y^{p-1}≡1(
yp−1≡1(mod
p
)
p)
p),所以暴力枚举只要枚举到
p
−
1
p-1
p−1即可。
但是这样还是不够优化,咋办?
令
x
=
m
i
−
j
x=mi-j
x=mi−j(其中
m
=
c
e
i
l
(
s
q
r
t
(
p
)
)
m=ceil(sqrt(p))
m=ceil(sqrt(p)),
c
e
i
l
ceil
ceil是向上取整的意思)所以呢,问题就变成了
y
m
i
−
j
≡
z
(
y^{mi-j}≡z(
ymi−j≡z(mod
p
)
p)
p),然后呢,问题就变成了
y
m
i
=
y
j
∗
z
(
y^{mi}=y^j*z(
ymi=yj∗z(mod
p
)
p)
p),好了,现在就好多了
接下来很容易发现
j
j
j<
m
m
m,所以我们可以暴力枚举
j
j
j,与所得
y
j
∗
z
(
y^j*z(
yj∗z(mod
p
)
p)
p)的存在哈希表里,然后再暴力枚举
m
i
mi
mi,最后得出结果
const int maxn=1000007;
int h[maxn],ne[maxn],bj[maxn],tot;
ll has[maxn];
void add(ll x,ll j) {
int kl=x%maxn;
tot++;bj[tot]=j;has[tot]=x;
ne[tot]=h[kl];h[kl]=tot;
}
ll find(ll x) {
ll kl=x%maxn,i;
for(i=h[kl];i!=-1;i=ne[i])
if(has[i]==x)return bj[i];
return -1;
}
ll bsgs(ll a,ll b,ll n) {//a^x=b(mod n)
if(b==1)return 0;
memset(h,-1,sizeof(h));
ll m=ceil(sqrt((double)n)),q=1,j,x=1;
for(ll i=0;i<m;i++) {add((q*b)%n,i);q=(q*a)%n;}
for(ll i=m;;i+=m) {
x=(x*q)%n;j=find(x);
if(j!=-1) return i-j;
if(i>n) break;
}
return -1;
}
#扩展bsgs
扩展bsgs更加神奇和美妙~~因为——
扩展bsgs不要求p是质数,多美妙啊!
这是怎么做到的呢?
大家意会一下就知道,若
t
=
g
c
d
(
a
,
p
)
t=gcd(a,p)
t=gcd(a,p).
那么
y
≡
z
(
y≡z(
y≡z(mod
p
)
p)
p)等价于
a
t
≡
b
∗
t
(
at≡b*t(
at≡b∗t(mod
c
∗
t
)
c*t)
c∗t)
即
a
≡
b
(
a≡b(
a≡b(mod
c
)
c)
c)
得到下面这段代码(以下a表示原题中y,b表示原题中z):
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;
}
解说:
每次循环使得
a
x
−
1
∗
a
d
≡
b
d
(
a^{x-1}* \frac{a}{d} ≡ \frac{b}{d}(
ax−1∗da≡db(mod
p
d
)
\frac{p}{d})
dp)
所以最后
d
∗
a
x
−
c
n
t
≡
b
(
d*a^{x-cnt}≡b(
d∗ax−cnt≡b(mod
c
)
c)
c)(都是运行完循环后的)
那么最后只要设
x
−
c
n
t
=
i
∗
m
−
j
x-cnt=i*m-j
x−cnt=i∗m−j即可用普通bsgs求解
以我们做一下poj3243(QAQ哈希写错调半天)
代码:
#include<iostream>
#include<cstdio>
#include<climits>
#include<algorithm>
#include<cstring>
#include<cmath>
using namespace std;
#define ll long long
#define mod 99991
int h[mod],tot;
ll ne[mod<<1],bj[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 kl=num%mod;
bj[++tot]=x;has[tot]=num;ne[tot]=h[kl];h[kl]=tot;
}
ll find(ll num){
ll kl=num%mod;
for(ll i=h[kl];i!=-1;i=ne[i])
if(has[i]==num)return bj[i];
return -1;
}
ll solve(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;
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;
}
int main()
{
ll a,b,p,ans;
while(1){
scanf("%lld%lld%lld",&a,&p,&b);
if(!a&&!b&&!p)break;
memset(h,-1,sizeof(h));tot=0;
ans=solve(a,b,p);
if(ans==-1)puts("No Solution");
else printf("%lld\n",ans);
}
return 0;
}