再观BZOJ2219数论之神...

shinbokuow posted @ 10 年前 in BZOJ with tags 数论数学 BSGS 原根 EXBSGS , 1987 阅读

题解:

当时A掉的时候其实并不会...

于是现在膜别人的题解终于会啦!

求方程

$x^a\equiv{b}\pmod{n}$ 在

$[0,n-1]$ 内的解的个数.

我们考虑通过对于两边取指标,转化为线性同余方程问题.

但是发现

$n$ 不一定存在原根,我们将n进行质因数分解.

原方程的解的个数等于每一个子方程

$x^a\equiv{b}\pmod{p_i^{q_i}}$ 解的个数的乘积.

我们考虑中国剩余定理,就能发现这样的乘法原理显然成立.

于是我们只考虑每一个子方程的解的个数.

我们考虑以下几种情况:

(1)

$b$ 存在指标,即

$(b,p_i)=1$

我们对两边取指标得到

$a\times{ind(x)}\equiv{ind(b)}\pmod{\phi(p_i^{q_i})}$

由线性同余方程的知识可知,如果

$(a,\phi(p_i^{q_i}))\not|ind(b)$ ,则无解;

否则解的个数为

$(a,\phi(p_i^{q_i}))$ .

(2)

$b$ 不存在指标,且

$b=0$

我们可以得到

$p_i|x$ ,不妨令

$x=cp_i^{x'}((c,p_i)=1)$ ,可以得到

$ax'\geq{q_i},即x'\geq{\lceil\frac{q_i}{a}\rceil}$ .

所以我们有只需满足

$p_i^{\lceil\frac{q_i}{a}\rceil}|x$ ,于是解的个数就为

$p_i^{q_i-\lceil\frac{q_i}{a}\rceil}$ .

(3)

$b$ 不存在指标,且

$b\not=0$

不妨令

$b=cp_i^{b'}((c,p_i)=1)$ .

若

$b'\%a\not=0$ ,显然无解.

否则我们将等式两边同除

$p_i^{b'}$ 得到:

$(\frac{x}{p_i^{\frac{b'}{a}}})^a\equiv{c}\pmod{p_i^{q_i-b'}}$

令

$x'=\frac{x}{p_i^{\frac{b'}{a}}}$ ,由于满足(2)的条件,我们可以求出

$x'$ 的解的个数.

注意

$x'$ 的范围是

$[0,p_i^{q_i-b'})$ ,但原式的

$\frac{x}{p_i^{\frac{b'}{a}}}$ 范围是

$[0,p_i^{q_i-\frac{b'}{a}})$ .

这相当于我们求出的解是在每个长度为

$p_i^{q_i-b'}$ 的区间中的解的个数.

所以我们需要将解的个数乘以

$p^{b'-\frac{b'}{a}}$ .

于是我们就彻底解决了这个问题!

代码:

100

101

102

103

104

105

106

107

108

109

110

111

112

113

114

115

116

117

118

119

120

121

122

123

124

125

126

127

128

129

130

131

132

133

134

135

136

137

138

139

140

141

142

143

144

145

146

147

148

149

150

151

152

153

154

155

156

157

158

159

160

161

162

163

164

165

166

167

168

169

170

171

172

173

174

175

176

177

#include<cmath>
#include<cstdio>
#include<cctype>
#include<climits>
#include<cstring>
#include<cstdlib>
#include<iostream>
#include<algorithm>
using namespace std;
  
typedef long long ll;
  
ll pow(ll x,ll y,ll p){
    ll re=1,t=x;
    for(;y;y>>=1,t=t*t%p)
        if(y&1)
            re=re*t%p;
    return re;
}
ll pow(ll x,ll y){
    ll re=1,t=x;
    for(;y;y>>=1,t*=t)
        if(y&1)
            re*=t;
    return re;
}
  
ll gcd(ll a,ll b){
    return(!b)?a:gcd(b,a%b);
}
  
inline void exgcd(ll a,ll b,ll&d,ll&x,ll&y){
    if(!b){
        d=a;
        x=1;
        y=0;
    }
    else{
        exgcd(b,a%b,d,y,x);
        y-=x*(a/b);
    }
}
inline ll inv(ll a,ll n){
    ll d,x,y;
    exgcd(a,n,d,x,y);
    return(x%n+n)%n;
}
  
inline ll get_g(ll p){
    static int pr[110];
    int top=0;
    ll t=p-1;
    for(int i=2;i*i<=p-1;++i){
        if(t%i==0){
            pr[++top]=i;
            while(t%i==0)
                t/=i;
        }
    }
    if(t!=1)
        pr[++top]=t;
    for(int i=1;;++i){
        bool ok=1;
        for(int j=1;j<=top;++j)
            if(pow(i,(p-1)/pr[j],p)==1){
                ok=0;
                break;
            }
        if(ok)
            return i;
    }
}
  
struct Hashset{
    static const int mod=13131;
    int head[mod],next[30010],f[30010],v[30010],ind;
    inline void reset(){
        ind=0;
        memset(head,-1,sizeof head);
    }
    inline void insert(int x,int _v){
        int ins=x%mod;
        for(int j=head[ins];j!=-1;j=next[j])
            if(f[j]==x){
                v[j]=min(v[j],_v);
                return;
            }
        f[ind]=x;
        v[ind]=_v;
        next[ind]=head[ins];
        head[ins]=ind++;
    }
    inline int operator[](int x){
        int ins=x%mod;
        for(int j=head[ins];j!=-1;j=next[j])
            if(f[j]==x)
                return v[j];
        return -1;
    }
}M;
  
ll BSGS(ll c,ll a,ll b,ll p){
    ll m=ceil(sqrt(p)),t=1,ct=c;
    M.reset();
    for(int i=0;i<m;++i,(t*=a)%=p,(ct*=a)%=p)
        M.insert(ct,i);
    t=inv(t,p);
    for(int i=0;i*m<p;++i,(b*=t)%=p)
        if(M[b]!=-1)
            return i*m+M[b];
    return -1;
}
ll EXBSGS(ll a,ll b,ll p){
    ll c=1,tim=0,_gcd;
    while((_gcd=gcd(a,p))!=1){
        if(b%_gcd)
            return -1;
        p/=_gcd;
        ++tim;
        (c*=(a/_gcd))%=p;
        (b/=_gcd)%=p;
    }
    ll ans=BSGS(c,a,b,p);
    if(ans!=-1)
        ans+=tim;
    return ans;
}
  
ll solve(ll a,ll b,ll p,ll t){
    b%=pow(p,t);
    if(b==0)
        return pow(p,t-((t-1)/a+1));
    int _b=0;
    while(b%p==0)
        b/=p,++_b;
    ll phi=pow(p,t)*(p-1)/p;
    ll g=get_g(p);
    if(_b==0){
        ll ind_b=EXBSGS(g,b,pow(p,t));
        if(ind_b%gcd(a,phi))
            return 0;
        return gcd(a,phi);
    }
    else{
        if(_b%a)
            return 0;
        else
            return solve(a,b,p,t-_b)*pow(p,_b-_b/a);
    }
}
  
int main(){
#ifndef ONLINE_JUDGE
    freopen("tt.in","r",stdin);
#endif
    int T;
    cin>>T;
    ll A,B,K,n,_n;
    int cnt;
    while(T--){
        scanf("%lld%lld%lld",&A,&B,&K);
        n=_n=K<<1^1;
        ll ans=1;
        for(int i=2;i*i<=n;++i){
            if(_n%i==0){
                cnt=0;
                while(_n%i==0)
                    ++cnt,_n/=i;
                ans*=solve(A,B,i,cnt);
            }
        }
        if(_n!=1)
            ans*=solve(A,B,_n,1);
        printf("%lld\n",ans);
    }
    return 0;
}

[回复]

ww140142 说:
10 年前

看个博客还有注册个号差评。。

[回复]

wyfcyx 说:
10 年前

@ww140142: 那是为了要评论没办法啊。。。

[回复]

ww140142 说:
10 年前

并不是这个意思，我是说那个网址进去要在OJ注册账号= =

[回复]

wyfcyx 说:
10 年前

@ww140142: 好吧0.0

再观BZOJ2219数论之神...

search

lastest comments

friend links

visitor counter