# 【BZOJ3512】DZY Loves Math IV（杜教筛）

### 推式子

$S_n(m)=\sum_{i=1}^m\phi(ni)$

$S_n(m)=p\times \sum_{i=1}^m\phi(qi)=p\times \sum_{i=1}^m\phi(\frac q{gcd(q,i)})\phi(i)\times gcd(q,i)$

$S_n(m)=p\times \sum_{i=1}^m\phi(\frac q{gcd(q,i)})\phi(i)\times\sum_{d|gcd(q,i)}\phi(d)$

$S_n(m)=p\times\sum_{i=1}^m\phi(i)\sum_{d|gcd(q,i)}\phi(\frac qd)$

$S_n(m)=p\times \sum_{d|q}\phi(\frac qd)\sum_{i=1}^{\lfloor\frac md\rfloor}\phi(di)=p\times \sum_{d|q}\phi(\frac qd)S_d(\lfloor\frac md\rfloor)$

### 代码

#pragma GCC optimize(2)
#pragma GCC optimize("inline")
#include<bits/stdc++.h>
#define Tp template<typename Ty>
#define Ts template<typename Ty,typename... Ar>
#define Reg register
#define RI Reg int
#define Con const
#define CI Con int&
#define I inline
#define W while
#define N 100000
#define X 1000000007
#define Inc(x,y) ((x+=(y))>=X&&(x-=X))
using namespace std;
int n,m,sm;
class DuSieve//杜教筛
{
private:
#define SZ 100000
#define V(x) ((x)<=sm?f1[x]:f2[m/(x)])
int Pt,P[SZ+5],phi[SZ+5],sphi[N+5],f1[SZ+5],f2[SZ+5];map<int,int> s[SZ+5];
I int SPhi(CI x)//求Phi
{
if(x<=SZ) return sphi[x];if(V(x)) return V(x);RI l,r,t=1LL*x*(x+1)/2%X;
for(l=2;l<=x;l=r+1) r=x/(x/l),t=(t-1LL*(r-l+1)*SPhi(x/l)%X+X)%X;return V(x)=t;
}
public:
I DuSieve()//预处理
{
phi[1]=1;for(RI i=2,j,t;i<=SZ;++i)
for(!P[i]&&(phi[P[++Pt]=i]=i-1),j=1;j<=Pt&&(t=i*P[j])<=SZ;++j)
if(P[t]=1,i%P[j]) phi[t]=phi[i]*(P[j]-1);else {phi[t]=phi[i]*P[j];break;}
for(RI i=1;i<=SZ;++i) sphi[i]=(sphi[i-1]+phi[i])%X;
}
I int S(RI n,CI m)//求S
{
if(!m) return 0;if(n==1) return SPhi(m);if(m==1) return phi[n];//边界
if(s[n][m]) return s[n][m];RI i,k=n,p=1,q=1,t=0;//记忆化
for(i=1;P[i]*P[i]<=k;++i) if(!(k%P[i])) for(k/=P[i],q*=P[i];!(k%P[i]);k/=P[i],p*=P[i]);//质因数分解
for(k^1&&(q*=k),i=1;i*i<=q;++i) !(q%i)&&//暴力枚举q的因数
(t=(1LL*phi[q/i]*S(i,m/i)+t)%X,i^(q/i)&&(t=(1LL*phi[i]*S(q/i,m/(q/i))+t)%X));
return s[n][m]=1LL*p*t%X;//注意最后乘上p
}
}D;
int main()
{
RI i,j,t=0;for(scanf("%d%d",&n,&m),sm=sqrt(m),i=1;i<=n;++i) Inc(t,D.S(i,m));//暴力枚举n
return printf("%d\n",t),0;
}

posted @ 2020-06-15 20:40  TheLostWeak  阅读(15)  评论(0编辑  收藏