zoukankan      html  css  js  c++  java
  • 数论细节梳理&模板

    初阶

    扩展欧拉

    (kgevarphi(m))时,(b^kequiv b^{k\%varphi(m)+varphi(m)}(mod m))

    扩展CRT

    推式子合并同余方程。
    https://www.cnblogs.com/flashhu/p/9346078.html

    扩展BSGS

    根号大暴力,就是细节很多。
    https://www.cnblogs.com/flashhu/p/9737769.html

    扩展Lucas

    对于模数不是质数的,唯一分解为(prod p_i^{k_i}),算出组合数对每一个(p_i^{k_i})取模的结果,用CRT合并。
    问题在于求(n!\%p^k)
    (n!)中所有(p)因子提出来,因子个数(sumlimits_{i=1}lfloorfrac{n}{p^i} floor)
    剩下的数的贡献单独考虑,记为(F_n),发现它在模意义下以(p^k)为周期,(fac_n=F_{p^k}^{lfloorfrac{n}{p^k} floor}F_{n\%p^k}fac_{lfloorfrac n p floor}),递归处理。
    时间复杂度(sum p_i^{k_i})
    洛谷P4720 【模板】扩展卢卡斯

    #include<bits/stdc++.h>
    #define LL long long
    #define R register LL
    using namespace std;
    LL n,m,YL,F[1000009];
    void exgcd(R a,R b,R&x,R&y){
        if(!b){x=1;y=0;return;}
        exgcd(b,a%b,y,x);y-=a/b*x;
    }
    inline LL qpow(R b,R k,R p){
        R a=1;
        for(;k;k>>=1,b=b*b%p)
            if(k&1ll)a=a*b%p;
        return a;
    }
    LL fac(R n,R p,R k){
        return n?qpow(F[k],n/k,k)*F[n%k]%k*fac(n/p,p,k)%k:1;
    }
    inline LL inv(R n,R p){
        R x,y;exgcd(n,p,x,y);
        return x<=0?x+p:x;
    }
    inline LL cnt(R n,R p){
        R k=0;
        for(n/=p;n;n/=p)k+=n;
        return k;
    }
    inline LL C(R n,R m,R p,R k){
        for(R i=F[0]=1;i<=k;++i)
            F[i]=i%p?F[i-1]*i%YL:F[i-1];
        return fac(n,p,k)*inv(fac(m,p,k),k)%k*inv(fac(n-m,p,k),k)%k*qpow(p,cnt(n,p)-cnt(m,p)-cnt(n-m,p),k)%k;
    }
    inline LL CRT(R a,R p){
        return YL/p*inv(YL/p,p)%YL*a%YL;
    }
    int main(){
        cin>>n>>m>>YL;
        R now=YL,lim=sqrt(YL),ans=0;
        for(R p=2;p<=lim;++p){
            if(now%p)continue;
            R k=1;while(now%p==0)now/=p,k*=p;
            ans=(ans+CRT(C(n,m,p,k),k))%YL;
            lim=sqrt(now);
        }
        if(now>1)ans=(ans+CRT(C(n,m,now,now),now))%YL;
        cout<<ans<<endl;
        return 0;
    }
    

    进阶

    Miller-Rabin&Pollard-Rho

    背板子
    洛谷P4718 【模板】Pollard-Rho算法
    辣鸡卡常题,参考了一下讨论版里的卡常技巧
    upd:模数为long long级别的乘法取模原理
    如果NOI能用什么__int128的话倒还好说
    但是别忘了long double__int128的位长是一样的啊!因此用它做long long级别的运算的精度损失是完全可以接受的。
    我们知道(a*bmod p=a*b-lfloorfrac{a*b}p floor*p)
    然后我们用long double先除后乘,考虑浮点数除法会丢掉小数点后面一个很小的值,而之后强制转换整形是要强制舍去的,所以要加一个eps。

    #include<bits/stdc++.h>
    #define SL __int128
    #define LL long long
    #define RG register
    #define R RG int
    using namespace std;
    LL t,n,ans;
    inline LL Mul(LL a,LL b,LL p){//㧟的玄学乘法取模
        LL d=((long double)a/p*b+1e-8);
        LL r=a*b-d*p;
        return r<0?r+p:r;
    }
    inline LL Gcd(LL a,LL b){//㧟的builtin辗转相减
        if(!a||!b)return a|b;
        int t=__builtin_ctzll(a|b);
        a>>=__builtin_ctzll(a);
        do{
            b>>=__builtin_ctzll(b);
            if(a>b)swap(a,b);
            b-=a;
        }while(b);
        return a<<t;
    }
    inline LL Qpow(LL b,LL k,LL p,LL a=1){
        for(;k;k>>=1,b=(SL)b*b%p)
            if(k&1ll)a=(SL)a*b%p;
        return a;
    }
    inline LL MR(LL n){
        if(n==2)return 1;
        if(n==1||(1&n)==0)return 0;
        static int a[]={2,3,7,13,61,24251};
        LL u=n-1,x,y;R k=0;
        while((1&u)==0)u>>=1,++k;
        for(R i=0;i<6&&a[i]<n;++i){
            x=Qpow(a[i],u,n);
            for(R j=0;j<k;++j,x=y){
                y=(SL)x*x%n;
                if(y==1&&x!=1&&x!=n-1)return 0;
            }
            if(x!=1)return 0;
        }
        return 1;
    }
    inline LL F(LL a,LL c,LL n){
        LL t=Mul(a,a,n)+c;
        return t<n?t:t-n;
    }
    inline LL PR(LL n){
        if((1&n)==0)return 2;
        LL c=rand(),a=rand(),b=a,g;
        do{
            a=F(a,c,n);b=F(F(b,c,n),c,n);
            g=Gcd(n,abs(a-b));
            if(g!=1&&g!=n)return g;
        }while(a!=b);
        return g;
    }
    void find(LL n){
        if(n==1||n<=ans)return;
        if(MR(n)){ans=max(ans,n);return;}
        LL d=n;
        while(d==n)d=PR(n);
        while(n%d==0)n/=d;
        find(n);find(d);
    }
    int main(){
        srand(time(NULL));
        cin>>t;
        while(t--){
            cin>>n;ans=1;find(n);
            if(ans==n)puts("Prime");
            else cout<<ans<<'
    ';
        }
        return 0;
    }
    

    狄利克雷卷积&杜教筛

    积性函数的性质:两个积性函数的狄利克雷卷积还是积性函数。
    有一些推式子题里需要通过构造积性函数来加快函数求值。
    常用的卷积式
    (mu* extbf1=epsilonLeftrightarrowsumlimits_{d|n}mu(d)=[n==1])
    (varphi* extbf1= extbf{id}Leftrightarrowsumlimits_{d|n}varphi(d)=n)
    (mu* extbf{id}=varphiLeftrightarrowsumlimits_{d|n}frac{n}{d}mu(d)=varphi(n))
    ( extbf1* extbf1= extbf{d}; extbf{d}* extbf1=sigma)(感觉还有很多都是可以相互推出来的)

    杜教筛:求积性函数(f)的前缀和(记为(S))。
    构造积性函数(g),前提是(g)(F=f*g)的前缀和都可以(O(1))算。
    (sumlimits_{i=1}^nF(i)=sumlimits_{i=1}^nsumlimits_{xy=i}f(x)g(y)=sumlimits_{y=1}^ng(y)sumlimits_{x=1}^{lfloorfrac{n}{y} floor}f(x)=g(1)S(x)+sumlimits_{y=2}^ng(y)S(lfloorfrac n y floor))
    (g(1)S(x)=sumlimits_{i=1}^nF(i)-sumlimits_{y=2}^ng(y)S(lfloorfrac n y floor))
    看到了可以整除分块+递归的地方。
    预处理(F,g)(n^frac 2 3)的前缀和,时间复杂度(O(n^frac 2 3))
    洛谷P4213 【模板】杜教筛(Sum)

    #include<bits/stdc++.h>
    #define LL long long
    #define UI unsigned int
    #define RG register
    #define R RG int
    using namespace std;
    const LL N=2147483648,M=1.7e6,L=2000;
    int pr[M],cnt;bool vis[M];
    struct Dat{LL p;int u;}s[M],t[L];
    Dat&S(UI n){
        if(n<M)return s[n];
        UI x=N/n;if(vis[x])return t[x];vis[x]=1;
        Dat&ans=t[x];
        ans.p=(LL)n*(n+1)>>1;ans.u=1;
        for(UI l=2,r;l<=n;l=r+1){
            r=n/(n/l);Dat&ret=S(n/l);
            ans.p-=(r-l+1)*ret.p;
            ans.u-=(r-l+1)*ret.u;
        }
        return ans;
    }
    int main(){
        s[1]=(Dat){1,1};vis[1]=1;
        for(R i=2;i<M;++i){
            if(!vis[i])s[pr[++cnt]=i]=(Dat){i-1,-1};
            for(R j=1,x;j<=cnt&&(x=i*pr[j])<M;++j){
                vis[x]=1;
                if(i%pr[j]==0){s[x].p=s[i].p*pr[j];break;}
                s[x].p=s[i].p*(pr[j]-1);s[x].u=-s[i].u;
            }
            s[i].p+=s[i-1].p;s[i].u+=s[i-1].u;
        }
        memset(vis,0,M);
        R t,n;cin>>t;
        while(t--){
            cin>>n;
            Dat ret=S(n);
            cout<<ret.p<<' '<<ret.u<<endl;
            memset(vis,0,L);
        }
        return 0;
    }
    

    Min_25筛

    同样是求积性函数前缀和,不过扩展性很强。
    yyb巨佬告诉我这是今年新鲜出炉的筛法,吊打了某洲阁筛(话说写这句话的时候是今年的最后一天了)
    核心思想是体现在Step2中的将最小质因子相同的数的贡献一起算。
    应用前提:对于(pin P)(f(p))是一个简单多项式,(f(p^e))可以快速计算。比某杜教筛不知道高到哪里去了。
    时间复杂度玄学,但是跑得快
    https://www.cnblogs.com/cjyyb/p/10169190.html
    https://www.cnblogs.com/cjoieryl/p/9403579.html
    https://www.cnblogs.com/zhoushuyu/p/9187319.html
    https://www.cnblogs.com/cx233666/p/10173977.html
    https://www.cnblogs.com/GuessYCB/p/10061411.html (都是高级操作)

    Step1

    把多项式拆开,对每个(f(x)=x^k)的项分开算贡献。
    (sqrt N)范围的质数筛出来,记为(P)。预处理质数的函数前缀和(fs)
    先对于每个(ile sqrt N)(sumlimits_{x=1}^{lfloorfrac N i floor}[xin P]f(x))
    (g(n,j))(xin[1,n])中满足(xin P)(x)的最小质因子小于(P_j)(f(x))之和
    于是显然有(g(n,|P|)=sumlimits_{x=1}^{n}[xin P]f(x))
    实现过程:滚动数组,外层枚举(j),直接用一维数组保存每一个(n=lfloorfrac N i floor)当前的(g(n,j))
    内层枚举(i),有(g(n,j)=g(n,j-1)-f(P_j)cdot(g(lfloorfrac{n}{P_j} floor,j-1)-fs(j-1))),结合定义用容斥思想理解。

    Step2

    开始算答案了。
    (S(n,j))(xin[1,n])中满足(x)的最小质因子大于等于(P_j)(f(x))之和(注意和(g)的区别)
    于是我们要求的答案即是(S(N,1)+f(1))(S)里面不包含(f(1))
    (S(n,j))时,分质数和合数算贡献
    质数:(g(n,|P|)-fs(j-1))
    合数:枚举最小质因子(P_j)以及它的次数(e),相同的一起算(注意不包含(f(1)),因此还有一个尾项)
    (sumlimits_{k=j}^{P_k^2le n}sumlimits_{e=1}^{P_k^{e+1}le n}S(lfloorfrac{n}{P_k^e} floor,k+1)f(P_k^e)+f(P_k^{e+1}))
    LOJ6053 简单的函数

    #include<bits/stdc++.h>
    #define LL long long
    #define RG register
    #define R RG LL
    using namespace std;
    const LL N=2e5+9,YL=1e9+7;
    LL n,Sq,p,m,pr[N],id1[N],id2[N],w[N],fx[N],gx[N],g1[N];
    bool np[N];
    inline LL Getid(R x){
    	return x<=Sq?id1[x]:id2[n/x];
    }
    void Sieve(){
    	for(R i=2;i<=Sq;++i){
    		if(np[i])continue;
    		pr[++p]=i;fx[p]=(fx[p-1]+i)%YL;
    		for(R j=i*i;j<=Sq;j+=i)np[j]=1;
    	}
    }
    LL S(R n,R j){
    	if(n<=1||pr[j]>n)return 0;
    	R x=Getid(n);
    	R ret=(gx[x]-fx[j-1]-g1[x]+j-1+2*YL)%YL;
    	for(R k=j;k<=p&&pr[k]*pr[k]<=n;++k){
    		R p1=pr[k],p2=p1*p1;
    		for(R e=1;p2<=n;++e,p1=p2,p2*=pr[k])
    			ret=(ret+S(n/p1,k+1)*(pr[k]^e)+(pr[k]^(e+1)))%YL;
    	}
    	return ret;
    }
    int main(){
    	cin>>n;
    	if(n==1)return puts("1"),0;
    	Sq=sqrt(n);Sieve();
    	for(R i=1,j;i<=n;i=j+1){
    		R x=w[++m]=n/i;j=n/x;
    		(x<=Sq?id1[x]:id2[n/x])=m;
    		gx[m]=(x+2)%YL*((x-1)%YL)%YL*((YL+1)>>1)%YL;
    		g1[m]=(x-1)%YL;
    	}
    	for(R j=1;j<=p;++j)
    		for(R i=1;i<=m&&pr[j]*pr[j]<=w[i];++i){
    			R x=Getid(w[i]/pr[j]);
    			gx[i]=((gx[i]-pr[j]*(gx[x]-fx[j-1]))%YL+YL)%YL;
    			g1[i]=(g1[i]-g1[x]+j-1+YL)%YL;
    		}
    	return cout<<S(n,1)+3<<endl,0;
    }
    
  • 相关阅读:
    Django 之 ORM
    python 变量问题 提问有大佬知道吗?
    celery 使用详解
    ubuntu 网络图标消失问题解决
    django + uwsgi 部署上线
    Python pip常用源地址
    ConfigParser 模块 使用
    Mysql 远程登录问题解决
    商城项目 从设想到实现
    硬盘上有打开的锁和感叹号标志如何解决(win10系统)
  • 原文地址:https://www.cnblogs.com/flashhu/p/10190408.html
Copyright © 2011-2022 走看看