BZOJ 2219 数论之神 BSGS+CRT

题意:链接

方法: BSGS+CRT

解析:

这道题有什么区别呢?

就是他取模的值不是一个质数了

这怎么办呢?

我们来把这个数分解质因数。

P=∏piti

然后对于每一个piti我们单独计算。

最后用中国剩余定理合并

这其实跟想象中的中国剩余定理没有什么关系,其实是其中的一个性质,就是mod p的解数总数量等于每一个mod piti的乘积。

这是一个性质,别问我证明。

原方程xA=B(mod C)

转化为xA=B(mod Ci)?>Ci是质数或是质数的幂次。

我们可以单独计算。

这跟之前的题类似。

不过有一点需要注意的地方。

就是如果直接跑最后求乘积的话是错的。

为什么呢,因为这个解的范围显然改变了。原来为C,不过现在变成了每个Ci的范围,显然不对。

况且,对于B来说,其可能比Ci要大,显然是错误的。

那么对于一个Ci,我们需要知道什么呢,知道它的质因子Pi,幂次Ti

首先先将这个B取余一下Ci,让这个方程正常一点。

然后这里就有问题了。

如果B整除Ci怎么办呢?

整除Ci的话,那么我们不妨设所有答案的gcd为Pik,明确一下,如果要是解最多的话显然我们让这个k越小越好。

回归下方程xA=0(mod PiTi)

PiA?k=0(mod PiTi)

A?k>=Ti

所以这个k最小是多少呢?Ti?1A+1

这应该是比较好想的?

然后这种情况我们就讨论完了。

现在来讨论B不为零的时候。

这时候我们可以把B换一个形式。

比如Pip?W

然后呢,此时等式变成了什么

xA=Pip?W(mod PiTi)

既然这样的话,右边的Pip一定可以均分到左边对吧?

如果不均分的话,显然二者是不会同余的。

所以不均分的话显然该方程无解。

然后我们来考虑均分的情况

对于该方程,我们可以把左右及模数共有的部分给约去。

所以变成了yA=W(mod PiTi?k)且gcd(y,Pi)=1

显然所有的子方程都可以化成这个东西。

然后来考虑这个东西怎么求解。

设mod=PiTi?k

然后我们可以求出模数的原根g;

然后呢我们可以求出来W在模mod下对g的指标。

不妨设为ind_W;

我们再设y在模mod下对g的指标为indy

则原方程转化为

A?indy=indW(mod φ(mod))

为什么取原根的优越性就体现出来吧?

因为该方程的取余后的值是可以覆盖[1,φ(mod))的所有值,也就是可以将这个方程进一步转化。

至于这个方程怎么解?

很简单吧,这就是Ax=B(mod C)嘛,exgcd秒解。

但是这里有一个问题,也就是说我们如此解出来的方案的数量是否与原方程相同呢?

xA=Pip?W(mod PiTi)

yA=W(mod PiTi?k)且gcd(y,Pi)=1

就这两个方程,我们看到,差异性就体现在模的数部分,这会有什么影响呢?

别忘了,如果我们求出来所有的解的话,正常的话,应该是在(mod PiTi)意义下的,但是现在变成了(mod PiTi?k)意义下的解..所以我们显然会少计算解,少计算了多少呢?

先看我们把x缩小了多少

显然是y显然是xPipA这么多,所以我们解出来的东西实际应该在(mod PiTi?pA)这个范围内但是我们解的范围是多少呢?PiTi?k这个范围内的,并且我们知道,k=p即我们求出了PiTi?p这个范围内的

所以我们要将该方程最终解中乘以个PiTi?pAPiTi?p即Pip?pA。

然后我们解完所有的方程后合并就是最终解。

这题解我写的好累啊。

popoqqq爷有个例子,可以去看看,但是有点bug,细心的人自己发现吧。

代码:

#include <cmath>
#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#define N 140345
#define INF 0x7fffffffffffffffll
using namespace std;
typedef long long ll;
int T,cnt,tot,totc;
int head[N+10];
ll a,b,c;
ll prime[N];
ll the_prime_c[N];
ll cnt_num_c[N];
struct node
{
    ll from,to,val;
    int next;
}edge[N+10];
void cutit(ll x)
{
    totc=0;
    for(ll i=2;i*i<=x;i++)
    {
        if(x%i==0)
        {
            the_prime_c[++totc]=i;
            cnt_num_c[totc]=1;
            x/=i;
            while(x%i==0){x/=i;cnt_num_c[totc]++;}
        }
    }
    if(x!=1)the_prime_c[++totc]=x,cnt_num_c[totc]=1;
}
void init()
{
    memset(head,-1,sizeof(head));
    cnt=1;
}
void edgeadd(ll from,ll to,ll val)
{
    edge[cnt].from=from,edge[cnt].to=to,edge[cnt].val=val;
    edge[cnt].next=head[from];
    head[from]=cnt++;
}
void exgcd(ll a,ll b,ll &x,ll &y,ll &gcd)
{
    if(!b)
    {
        x=1,y=0,gcd=a;
        return;
    }
    exgcd(b,a%b,y,x,gcd);
    y=y-a/b*x;
}
void get_factor(ll x)
{
    tot=0;
    for(int i=2;i*i<=x;i++)
    {
        if(x%i==0)
        {
            prime[++tot]=i;
            while(x%i==0)x/=i;
        }
    }
    if(x!=1)prime[++tot]=x;
}
ll quick_my(ll x,ll y,ll mod)
{
    ll ret=1;
    while(y)
    {
        if(y&1)ret=(ret*x)%mod;
        x=(x*x)%mod;
        y>>=1;
    }
    return ret;
}
ll get_inv(ll x,ll mod)
{
    ll X,Y,GCD;
    exgcd(x,mod,X,Y,GCD);
    return (X%mod+mod)%mod;
}
ll get_phi(ll x)
{
    ll ret=x,tmp=x;
    for(ll i=2;i*i<=tmp;i++)
    {
        if(tmp%i==0)
        {
            ret=ret*(i-1)/i;
            while(tmp%i==0){tmp/=i;}
        }
    }
    if(tmp!=1)ret=ret*(tmp-1)/tmp;
    return ret;
}
int check(ll x,ll MOD,ll PHI)
{
    for(int i=1;i<=tot;i++)
    {
        if(quick_my(x,PHI/prime[i],MOD)==1)return 0;
    }
    return 1;
}
ll find_primitive_root(ll x)
{
    ll tmp=get_phi(x);
    get_factor(tmp);
    for(int i=2;;i++)
    {
        if(check(i,x,tmp))return i;
    }
}
ll Baby_Step_Giant_Step(ll A,ll B,ll C)
{
    init();
    ll m=(int)ceil(sqrt(C));
    ll k=1;
    for(int i=0;i<m;i++)
    {
        int flag=1;
        for(int j=head[k%N];j!=-1;j=edge[j].next)
        {
            if(edge[j].val==k){flag=0;break;}
        }
        if(flag)edgeadd(k%N,i,k);
        k=k*A%C;
    }
    ll invk=get_inv(k,C);
    ll invD=1;
    for(int i=0;i<=m;i++)
    {
        ll tmpB=B*invD%C;
        for(int j=head[tmpB%N];j!=-1;j=edge[j].next)
        {
            if(edge[j].val==tmpB)return edge[j].to+i*m;
        }
        invD=invD*invk%C;
    }
    return -1;
}
ll work(ll a,ll b,ll pri,ll num)
{
    ll B=quick_my(pri,num,INF);b%=B;
    if(b==0)return quick_my(pri,num-((num-1)/a+1),INF);
    ll cnt_b=0;
    while(b%pri==0)
    {
        cnt_b++,b/=pri,num--;
    }
    if(cnt_b%a!=0)return 0;
    ll g=find_primitive_root(B);
    ll indB=Baby_Step_Giant_Step(g,b,B);
    ll X,Y,GCD;
    exgcd(a,get_phi(B),X,Y,GCD);
    // A*ind_g(A)=ind_g(B)(mod φ(C1));
    if(indB%GCD!=0)return 0;
    return GCD*quick_my(pri,cnt_b-cnt_b/a,INF);
}
int main()
{
    scanf("%d",&T);
    while(T--)
    {
        scanf("%lld%lld%lld",&a,&b,&c);
        c=c*2+1;
        cutit(c);
        ll ans=1;
        for(int i=1;i<=totc;i++)
        {
            ans*=work(a,b,the_prime_c[i],cnt_num_c[i]);
        }
        printf("%lld\n",ans);
    }
}

版权声明:本文为博主原创文章,未经博主允许不得转载。

时间: 2024-10-14 14:22:15

BZOJ 2219 数论之神 BSGS+CRT的相关文章

BZOJ 2219: 数论之神

题目:http://www.lydsy.com/JudgeOnline/problem.php?id=2219 N次剩余+CRT... 就是各种奇怪的分类讨论.. #include<cstring> #include<iostream> #include<cstdio> #include<map> #include<cmath> #include<algorithm> #define rep(i,l,r) for (int i=l;i

BZOJ 2219 数论之神 数论

题目大意:求在[0,p)范围内的解的个数 鏼爷的题解:http://jcvb.is-programmer.com/posts/42036 我只是来粘代码的QAQ 指标啥的原根啥的中国剩余定理啥的真的完全不知道QAQ #include <cmath> #include <cstdio> #include <cstring> #include <iostream> #include <algorithm> #define INF static_cas

【BZOJ】【2219】数论之神

中国剩余定理+原根+扩展欧几里得 题解:http://blog.csdn.net/regina8023/article/details/44863519 新技能get√: 1 LL Get_yuangen(LL p,LL phi){ 2 int c=0; 3 for(int i=2;i*i<=phi;i++) 4 if (phi%i==0) 5 f[++c]=i,f[++c]=phi/i; 6 for(int g=2;;g++){ 7 int j; 8 for(j=1;j<=c;j++) if

Bzoj2219 数论之神

Time Limit: 3 Sec  Memory Limit: 259 MBSubmit: 954  Solved: 268 Description 在ACM_DIY群中,有一位叫做“傻崽”的同学由于在数论方面造诣很高,被称为数轮之神!对于任何数论问题,他都能瞬间秒杀!一天他在群里面问了一个神题: 对于给定的3个非负整数 A,B,K 求出满足 (1) X^A = B(mod 2*K + 1) (2) X 在范围[0, 2K] 内的X的个数!自然数论之神是可以瞬间秒杀此题的,那么你呢? Inpu

BZOJ 3122 [Sdoi2013]随机数生成器 BSGS

题意:链接 方法: BSGS 解析: 首先他给出了你数列在mod p意义下的递推式. 所以我们可以求出来通项. Xn+1+k=a?(Xn+k) 所以b=(a?1)?k 则我们可以解出来k 那么这个数列的通项是什么呢? Xn=an?1?(X1+k)?k 题中给定Xn 求出n就行了. 所以只需要移项就好了. 这里有个问题,此时我们的通项公式是不包含首项的,所以需要特判首项,另外还有第一项以外为常数项的时候. 代码: #include <cmath> #include <cstdio>

BZOJ 3239 Discrete Logging(BSGS)

[题目链接] http://www.lydsy.com/JudgeOnline/problem.php?id=3239 [题目大意] 计算满足 Y^x ≡ Z ( mod P) 的最小非负整数 [题解] BSGS裸题. [代码] #include <cstdio> #include <cmath> #include <map> #include <algorithm> #include <tr1/unordered_map> using name

bzoj 3239: Discrete Logging【BSGS】

BSGS的板子题 此时 \( 0 \leq x \leq p-1 \) 设 \( m=\left \lceil \sqrt{p} \right \rceil ,x=i*m-j \)这里-的作用是避免逆元 于是可以把式子变形成这样:\( a^{im}\equiv ba^j(mod p) \) 枚举右边\( 0 \leq j <m \) ,用map或者hash以模数为下标来存每一个j 枚举左边\( 0 \leq i <m \) ,在map或者hash中查找对应的模数 #include<ios

BZOJ 1037 生日聚会(神DP)

这题的DP很难想,定义dp[i][j][a][b]表示用了i个男生,j个女生,任一连续的后缀区间内,男生比女生最多多a人,女生比男生最多多b人. 转移就是显然了. # include <cstdio> # include <cstring> # include <cstdlib> # include <iostream> # include <vector> # include <queue> # include <stack&

数论之神

[题目描述] 对于给定的3个非负整数A,B,K求出满足下列条件的X的个数. (1) X^A = B(mod 2*K + 1) (2) X ∈[0, 2K] [输入描述] 第一行有一个正整数T,表示接下来的数据的组数( T <= 1000) 之后对于每组数据,给出了3个整数A,B,K (1 <= A, B <= 10^9, 1 <= K <= 5 * 10^8). [输出描述] 输出一行,表示答案. [输入样例] 3 213 46290770 80175784 3 462907