码迷,mamicode.com
首页 > 其他好文 > 详细

bzoj 3512: DZY Loves Math IV

时间:2018-07-14 17:42:13      阅读:204      评论:0      收藏:0      [点我收藏+]

标签:namespace   bit   find   com   tin   end   gcd   des   定义   

Description

给定n,m,求
技术分享图片
模10^9+7的值。

Solution

\(S(n,m)\) 表示 \(\sum_{i=1}^{m}\phi(n*i)\)
\(Ans=\sum_{i=1}^{n}S(i,m)\)
\(S(n,m)=\sum_{i=1}^{m}\phi(n*i)\)

如果 \(\mu(n)!=0\)
则有 \(\sum_{i=1}^{m}\phi(\frac{n}{gcd(i,n)})*\phi(i)*gcd(i,n)\) (因为要保证除完\(gcd\)之后,两数不含相同的质因子,所以 \(\mu(n)!=0\))
\(\sum_{i=1}^{m}\phi(\frac{n}{gcd(i,n)})*\phi(i)*\sum_{d|i,d|n}\phi(d)\)
因为第一项和第三项是互质的 , 所以可以合并.
\(\sum_{i=1}^{m}\phi(i)*\sum_{d|n,d|i}\phi(\frac{n}{d})\)
\(\sum_{d=1}^{n}\phi(\frac{n}{d})*\sum_{i=1}^{\frac{m}{d}}\phi(d*i)\)
\(\sum_{d=1}^{n}\phi(\frac{n}{d})*S(d,\lfloor\frac{m}{d}\rfloor)\)
递归处理即可

如果 \(\mu(n)=0\)
我们直接提出 \(n\) 的多出的质因子之积 \(a\),使得 \(\mu(\frac{n}{a})!=0\)
那么 \(S(n,m)=\sum_{i=1}^{m}\phi(n*i)\) 中也可以提出 \(a\) 了,因为相同的质因子只会被算一次
根据定义式 \(\phi(n)=n*\Pi p_i\),所以 \(a\) 唯一的贡献就是使前面的 \(n\) 乘了个 \(a\)
\(a*S(n,m)=\sum_{i=1}^{m}\phi(\frac{n}{a}*i)=a*S(\frac{n}{a},m)\)
递归处理即可

边界条件 \(m=1?\) 时,结果为 \(\phi(n)?\), \(n=1?\) 时,跑一个杜教筛就行了

#include<bits/stdc++.h>
using namespace std;
const int N=1e5+10,mod=1e9+7;
int n,prime[N],num=0,m,phi[N],mu[N],pre[N],la[N],s[N];bool vis[N];
void priwork(){
    phi[1]=s[1]=1;
    for(int to,i=2;i<N;i++){
        if(!vis[i])prime[++num]=i,phi[i]=i-1,mu[i]=-1,pre[i]=i;
        for(int j=1;j<=num && i*prime[j]<N;j++){
            vis[to=i*prime[j]]=1;pre[to]=prime[j];
            if(i%prime[j])phi[to]=phi[i]*(prime[j]-1),mu[to]=-mu[i];
            else {phi[to]=phi[i]*prime[j];break;}
        }
        s[i]=(s[i-1]+phi[i])%mod;
    }
    for(int i=2;i<=n;i++){
        if(mu[i])continue;
        int last=0,x=i;la[i]=1;
        while(x>1){
            if(pre[x]==last)la[i]*=pre[x];
            last=pre[x];x/=pre[x];
        }
    }
}
map<int,int>S[N],T;
inline int calc(int n){
    if(n<N)return s[n];
    if(T.find(n)!=T.end())return T[n];
    int ret=(1ll*n*(n+1)>>1)%mod;
    for(int i=2,r;i<=n;i=r+1){
        r=n/(n/i);
        ret=(ret-1ll*calc(n/i)*(r-i+1))%mod;
    }
    if(ret<0)ret+=mod;
    return T[n]=ret;
}
inline int solve(int n,int m){
    if(m==1)return phi[n];
    if(n==1)return calc(m);
    if(S[n].find(m)!=S[n].end())return S[n][m];
    if(!mu[n])return 1ll*la[n]*solve(n/la[n],m)%mod;
    int ret=0,lim=min(m,(int)sqrt(n));
    for(int i=1;i<=lim;i++){
        if(n%i==0){
            if(i*i!=n)ret=(ret+1ll*phi[n/i]*solve(i,m/i)+1ll*phi[i]*solve(n/i,m/(n/i)))%mod;
            else ret=(ret+1ll*phi[n/i]*solve(i,m/i))%mod;
        }
    }
    return S[n][m]=ret;
}
int main(){
  freopen("pp.in","r",stdin);
  freopen("pp.out","w",stdout);
  cin>>n>>m;
  priwork();
  int ans=0;
  for(int i=1;i<=n;i++)ans=(ans+solve(i,m))%mod;
  printf("%d\n",ans);
  return 0;
}

bzoj 3512: DZY Loves Math IV

标签:namespace   bit   find   com   tin   end   gcd   des   定义   

原文地址:https://www.cnblogs.com/Yuzao/p/8546115.html

(0)
(0)
   
举报
评论 一句话评论(0
登录后才能评论!
© 2014 mamicode.com 版权所有  联系我们:gaon5@hotmail.com
迷上了代码!