BZOJ 3512 Dzy loves math IV

Link

Solution

n很小,不妨枚举。
对于一个固定的n
\sum_{i=1}^m\varphi(ni)。不妨设n为无平方因子数,那么

\begin{aligned}
S(n,m)=&\sum_{i=1}^m\varphi(ni)\\
=&\sum_{i=1}^m\varphi(i)\varphi(\frac n{\gcd(i,n)})\gcd(i,n)\\
=&\sum_{i=1}^m\varphi(i)\varphi(\frac n{\gcd(i,n)})\sum_{d|\gcd(n,i)}\varphi(d)\\
=&\sum_{i=1}^m\varphi(i)\sum_{d|\gcd(n,i)}\varphi(\frac {\gcd(i,n)}d)\varphi(\frac n{\gcd(i,n)})\\
\end{aligned}

而由于n为无平方因子数,\gcd(\frac {n}{\gcd(i,n)},\gcd(i,n))=1
因此上式可以化为

\begin{aligned}
&\sum_{i=1}^m\varphi(i)\sum_{d|\gcd(n,i)}\varphi(\frac {\gcd(i,n)}d)\varphi(\frac n{\gcd(i,n)})\\
=&\sum_{i=1}^m\varphi(i)\sum_{d|\gcd(n,i)}\varphi(\frac nd)\\
=&\sum_{d|n}\varphi(\frac nd)\sum_{i=1}^{[\frac md]}\varphi(id)\\
=&\sum_{d|n}\varphi(\frac nd)\sum_{i=1}^{[\frac md]}\varphi(id)\\
=&\sum_{d|n}\varphi(\frac nd)S(d,[\frac md])
\end{aligned}

边界情况S(1,m)就用杜教筛求\varphi的前缀和。
而S的第一维状态数为n,第二维状态数为O(\sqrt m),因为第二维会出现的状态都是m经过多次除法然后下取整的结果,只有O(\sqrt m)种。而每次计算枚举n的约数复杂度为O(\sqrt n),因此总复杂度就是O(n(\sqrt m+\sqrt n)+T),其中T为调用杜教筛求\varphi的耗时。

Code

const int P=1e9+7;

int Add(int x,int const &y) {
    return (x+=y)>=P?x-P:x;
}

int Minus(int x,int const &y) {
    return (x-=y)<0?x+P:x;
}

int Mul(long long x,int const &y) {
    return x*y%P;
}

const int N=1e6,XN=N+11;

int prime[XN],pcnt,phi[XN],dnum[XN],sphi[XN];
void Prep() {
    static bool notPrime[XN];
    phi[1]=dnum[1]=1;
    for(int i=2;i<=N;++i) {
        if(!notPrime[i]) {
            prime[++pcnt]=i;
            dnum[i]=i;
            phi[i]=i-1;
        }
        for(int j=1;j<=pcnt && i*prime[j]<=N;++j) {
            notPrime[i*prime[j]]=1;
            if(i%prime[j]==0) {
                dnum[i*prime[j]]=dnum[i];
                phi[i*prime[j]]=phi[i]*prime[j];
                break;
            } else {
                dnum[i*prime[j]]=dnum[i]*prime[j];
                phi[i*prime[j]]=phi[i]*phi[prime[j]];
            }
        }
    }
    for(int i=1;i<=N;++i)
        sphi[i]=Add(sphi[i-1],phi[i]);
}

int Sphi(int n) {
    if(n<=N)
        return sphi[n];
    else {
        static std::unordered_map<int,int> map;
        if(map.count(n))
            return map[n];
        int &res=map[n]=(long long)n*(n+1)/2%P;
        for(int l=2,r;l<=n;l=r+1) {
            r=n/(n/l);
            res=Minus(res,Mul(r-l+1,Sphi(n/l)));
        }
        return res;
    }
}

int Phi(int n) {
    if(n<=N)
        return phi[n];
    else {
        static std::unordered_map<int,int> map;
        if(map.count(n))
            return map[n];
        int &res=map[n]=n;
        for(int d=2,m=sqrt(n)+0.5;d<=m;++d)
            if(n%d==0) {
                do n/=d;
                while(n%d==0);
                (res/=d)*=(d-1);
            }
        if(n!=1)
            (res/=n)*=(n-1);
        return res;
    }
}

int S(int n,int m) {
    if(n==1)
        return Sphi(m);
    else {
        static std::map<std::pair<int,int>,int> map;
        std::pair<int,int> pid(n,m);
        if(map.count(pid))
            return map[pid];
        int &res=map[pid]=0;
        for(int d=1;d*d<=n;++d)
            if(n%d==0) {
                res=Add(res,Mul(Phi(n/d),S(d,m/d)));
                if(d*d!=n)
                    res=Add(res,Mul(Phi(d),S(n/d,m/(n/d))));
            }
        return res;
    }
}

int main() {
    Prep();
    int n,m,Ans=0;fin>>n>>m;
    for(int i=1;i<=n;++i)
        Ans=Add(Ans,Mul(i/dnum[i],S(dnum[i],m)));
    fout<<Ans<<'\n';
    return 0;
}

发表评论

电子邮件地址不会被公开。 必填项已用*标注