zoukankan      html  css  js  c++  java
  • 【luogu3768】简单的数学题【杜教筛】【欧拉函数】

    题目传送门
    题意
    (i=1nj=1nijgcd(i,j)) mod pn<=1010p为质数。
    题解
    之前一直想着用μ函数搞,结果死都只能搞到O(n)。后来,我fa♂现了bzoj3518点组计数,这道题用了欧拉函数的一个性质n=i|nφ(i),好像叫做欧拉反演?
    于是新日暮里的大门打开了。
    于是开始推式子。
    i=1nj=1nijgcd(i,j)
    =i=1nj=1nijk|gcd(i,j)φ(k)
    =k=1nk2φ(k)i=1nkj=1nkij
    我们令S(n)=i=1ni=n(n+1)2
    则原式
    =k=1nk2φ(k)S(nk)2
    于是我们只需要能够快速求出k=1nk2φ(k)就可以分块求和计算了。
    这不就是杜教筛干的活吗?
    i=1ni2j|iφ(j)
    =i=1ni3
    =S(n)2
    =j=1ni=1nji2j2φ(j)
    因为inj
    所以jni
    因此原式
    =i=1ni2j=1nij2φ(j)
    i=1带人
    i=1ni2φ(i)=S(n)2i=2ni2j=1nij2φ(j)
    然后就递归计算即可。细节详见代码。
    这道题就解决了!
    代码

    #include<cstdio>
    typedef long long ll;
    const ll N=5000005,M=1000017;
    ll mod,n,ans,p[N],phi[N];
    bool vis[N];
    struct hashmap{
        int cnt,head[M],nxt[M*10];
        ll to[M*10],v[M*10];
        bool count(ll f){
            ll x=f%M;
            for(int i=head[x];i;i=nxt[i]){
                if(to[i]==f){
                    return true;
                }
            }
            return false;
        }
        ll get(ll f){
            ll x=f%M;
            for(int i=head[x];i;i=nxt[i]){
                if(to[i]==f){
                    return v[i];
                }
            }
            return 0;
        }
        void set(ll f,ll val){
            ll x=f%M;
            to[++cnt]=f;
            nxt[cnt]=head[x];
            v[cnt]=val;
            head[x]=cnt;
        }
    }mp;
    inline ll sqr(ll x){
        x%=mod;
        return x*x%mod;
    }
    ll calc(ll x){
        if(!x){
            return 0;
        }
        x%=mod;
        return sqr((1+x)*x/2)%mod;
    }
    ll get(ll n){
        n%=mod;
        if(n%3!=1){
            return n*(n+1)/6%mod*(2*n+1)%mod;
        }else if(n&1){
            return (n+1)*(2*n+1)/6%mod*n%mod;
        }else{
            return n*(2*n+1)/6%mod*(n+1)%mod;
        }
    }
    ll solve(ll n){
        if(n<=5000000){
            return phi[n];
        }
        if(mp.count(n)){
            return mp.get(n);
        }
        ll res=calc(n);
        for(ll i=2,last;i<=n;i=last+1){
            last=n/(n/i);
            res-=solve(n/i)*(get(last)-get(i-1))%mod;
            res%=mod;
        }
        mp.set(n,res);
        return res;
    }
    int main(){
        scanf("%lld%lld",&mod,&n);
        phi[1]=1;
        for(int i=2;i<=5000000;i++){
            if(!vis[i]){
                p[++p[0]]=i;
                phi[i]=i-1;
            }
            for(int j=1;j<=p[0]&&i*p[j]<=5000000;j++){
                vis[i*p[j]]=true;
                if(i%p[j]){
                    phi[i*p[j]]=phi[i]*(p[j]-1);
                }else{
                    phi[i*p[j]]=phi[i]*p[j];
                }
            }
        }
        for(int i=1;i<=5000000;i++){
            phi[i]=phi[i]*i%mod*i%mod;
            phi[i]+=phi[i-1];
            phi[i]%=mod;
        }
        for(ll i=1,last;i<=n;i=last+1){
            last=n/(n/i);
            ans+=(solve(last)-solve(i-1))*calc(n/i)%mod;
            ans%=mod;
        }
        printf("%lld
    ",(ans+mod)%mod);
        return 0;
    }

    最后附上蒟蒻博主龙hou飞yan凤wu舞chi的公式草稿一张2333 鼠标写字真是太方便啦
    这里写图片描述

  • 相关阅读:
    Python实现客观赋权法
    Python实现熵值法确定权重
    正则化项L1和L2
    特征工程的归一化和标准化
    CentOS下Neo4j安装教程
    Window下Neo4j安装教程
    Window下JDK安装教程
    Git 命令
    Kubernetes 资源清单 常用字段,Pod 实例
    kubernetes 集群搭建
  • 原文地址:https://www.cnblogs.com/2016gdgzoi471/p/9476851.html
Copyright © 2011-2022 走看看