本来是想和 Min_25 一起写的，但是太长了（为了避免莫反那篇一样的惨案）所以就分开了w.
0 - 前置
莫比乌斯反演
详见 莫比乌斯反演学习笔记 ，这里摘录一些式子。
Dirichlet卷积：
[(f*g)(x)=sum_{d|x}f(d)g(frac{x}{d})
]
莫比乌斯反演：
[f=g*1=>g=f*mu\
f(x)=sum_{d|x}g(d)=>g(x)=sum_{d|x}f(d) imes mu(frac{x}{d}).
]
欧拉函数
性质：(sum_{d|n} varphi(d)=n)
表示成 Dirichlet卷积 的形式：(varphi*I=id) .
卷上 (mu) 得
[varphi*I=id\
=>varphi*I*mu=id*mu\
=>varphi * epsilon=id*mu
]
所以有
[varphi(n)=sum_{d|n}mu(d)cdotfrac{n}{d}
]
故
[frac{varphi(n)}{n}=sum_{d|n}frac{mu(d)}{d}
]
1 - 杜教筛
Method
用途：以低于线性的复杂度计算积性函数前缀和 ，即计算：(sum_{i=1}^n f(i)) .
推导：
构造两个积性函数 (h,g) ，使得 (h=f*g) .
记 (S(n)=sum_{i=1}^n f(i)) ，那么有
[egin{aligned}
sum_{i=1}^nh(i) &=sum_{i=1}^nsum_{d|i}g(d)cdot f(frac{i}{d})\\
&=sum_{d=1}^ng(d)cdot sum_{i=1}^{lfloor n/d
floor} f(i)\\
&=sum_{d=1}^ng(d)cdot S(lfloor frac{n}{d}
floor )\\
&=g(1)cdot S(n)+sum_{d=2}^{n}g(d)cdot S(lfloorfrac{n}{d}
floor)\\
g(1)S(n)&=sum_{i=1}^nh(i)-sum_{d=2}^ng(d)cdot S(lfloor frac{n}{d}
floor )
end{aligned}
]
所以关键就在于找一个容易得到的 (h) .求得之后对后面整除分块就好了，复杂度 (mathcal{O}(n^{frac 23})) .
关于复杂度 （@George1123）
直接递归求是 (Theta(n^{frac{3}{4}})) ：
[egin{split} T(n)=&sqrt n+sumlimits_{i=2}^{sqrt n}left(T(i)cdot T(lfloorfrac{n}{i} floor) ight)\ ge&sqrt n+sumlimits_{i=2}^{sqrt n}left(sqrt icdot sqrt{lfloorfrac{n}{i} floor)} ight)\ ge&sqrt n+sumlimits_{i=2}^{sqrt n}2sqrt{sqrt n}\ =&n^{frac{3}{4}}\ end{split} ]线性筛求出前 (n^{frac23}) 个，然后再杜教筛，是 (Theta(n^{frac{2}{3}})) .
[egin{split} T(n)=&m+sumlimits_{i=2}^{lfloorfrac nm floor}sqrt{lfloorfrac ni floor}\ =&m+frac{n}{sqrt m}\ ge&2n^{frac{2}{3}}(=: operatorname{while} m=n^{frac{2}{3}})\ end{split} ]
Example
莫比乌斯函数
求
[S(n)=sum_{i=1}^nmu(i).
]
由于 (mu*I=epsilon) ，所以 令 (I=>g,epsilon=>h) ,根据上面的式子
[g(1)S(n)=sum_{i=1}^nh(i)-sum_{d=2}^ng(d)cdot S(lfloor frac{n}{d}
floor )
]
有
[egin{aligned}
1 imes S(n)&=sum_{i=1}^nepsilon(i)-sum_{d=2}^nS(lfloor frac{n}{d}
floor )\\
S(n) &= 1-sum_{d=2}^nS(lfloor frac nd
floor)
end{aligned}
]
欧拉函数
求
[S(n)=sum_{i=1}^nvarphi(i)
]
由于 (varphi *I=ID) ，所以
[egin{aligned}
S(n) &=sum_{i=1}^n i-sum_{d=2}^n S(lfloor frac{n}{d}
floor)\\
&=frac{n(n+1)}{2}-sum_{d=2}^nS(lfloor frac{n}{d}
floor)
end{aligned}
]
？函数
求
[S(n)=sum_{i=1}^nicdot varphi(i)
]
考虑 Dirichlet卷积的形式：
[sum_{d|n}(dcdot varphi(d))cdot g(frac{n}{d})
]
令 (ID=>g) ，有
[sum_{d|n}(dcdot varphi(d))cdot frac{n}{d}=sum_{d|n}ncdot varphi(d)=nsum_{d|n}varphi(d)=n^2
]
所以有
[S(n)=sum_{i=1}^ni^2-sum_{d=2}^ndcdot S(lfloor frac{n}{d}
floor)
]
？？函数
求
[sum_{i=1}^nmu(i)i^2
]
令 (g(i)=i^2) ，有
[g*f(n)=sum_{d|n}mu(d)d^2cdotBig(frac{n}{d}Big)^2=sum_{d|n}mu(d)=[n=1]
]
所以
[S(n)=sum_{i=1}^nh(i)-sum_{d=2}^nd^2cdot S(lfloor frac{n}{d}
floor )=1-sum_{d=2^n}d^2cdot S(lfloor frac{n}{d}
floor)
]
Code
大常数选手打的 模板
//Author: RingweEH
const int N=5e5+10;
int pri[N],tot=0,n;
ll mu[N],phi[N];
bool is[N];
void Sieve()
{
mu[1]=phi[1]=1; is[1]=1;
for ( int i=2; i<=N-10; i++ )
{
if ( !is[i] ) { pri[++tot]=i,mu[i]=-1,phi[i]=i-1; }
for ( int j=1; j<=tot && (i*pri[j]<=(N-10)); j++ )
{
is[i*pri[j]]=1;
if ( i%pri[j]==0 ) { mu[i*pri[j]]=0,phi[i*pri[j]]=phi[i]*pri[j]; break; }
mu[i*pri[j]]=-mu[i]; phi[i*pri[j]]=phi[i]*phi[pri[j]];
}
}
for ( int i=2; i<=(N-10); i++ )
mu[i]+=mu[i-1],phi[i]+=phi[i-1];
}
map<int,ll> mp_mu,mp_phi;
ll Sieve_Phi( int n )
{
if ( n<=N-10 ) return phi[n];
if ( mp_phi[n] ) return mp_phi[n];
ll _res=0;
for ( int l=2,r; l<=n; l=r+1 )
{
r=n/(n/l); _res+=(r-l+1)*Sieve_Phi(n/l);
if ( r>=2147483647 ) break;
}
ll res=(ull)n*(n+1ll)/2ll-_res;
mp_phi[n]=res; return res;
}
ll Sieve_Mu( int n )
{
if ( n<=N-10 ) return mu[n];
if ( mp_mu[n] ) return mp_mu[n];
ll res=1;
for ( int l=2,r; l<=n; l=r+1 )
{
r=n/(n/l); res-=(r-l+1)*Sieve_Mu(n/l);
if ( r>=2147483647 ) break;
}
return mp_mu[n]=res;
}
int main()
{
Sieve();
int T=read();
while ( T-- )
{
int n=read();
printf("%lld %lld
",Sieve_Phi(n),Sieve_Mu(n) );
}
return 0;
}
下面是例题。
例题：毒瘤的数学题
求
[left(sumlimits_{i=1}^nsumlimits_{j=1}^nijgcd(i,j)
ight)mod p
]
(1le nle 10^{10},5 imes10^{8}le ple 1.1 imes 10^{9}) .
Solution
[egin{aligned}
sum_{i=1}^nsum_{j=1}^nijgcd(i,j) &=sum_{d=1}^n d^3sum_{i=1}^{lfloor n/d
floor}sum_{j=1}^{lfloor n/d
floor} ij[gcd(i,j)=1]\\
&=sum_{d=1}^nd^3sum_{i=1}^{lfloor n/d
floor}sum_{j=1}^{lfloor n/d
floor}ijsum_{k|d}mu(k)\\
令S(n)=frac{n*(n+1)}2,&=sum_{d=1}^nd^3sum_{k|d}mu(k)k^2S(lfloor frac{n}{dk}
floor)^2\\
&=sum_{d=1}^nd^3sum_{k=1}^{lfloor n/d
floor}mu(k)k^2S(lfloor frac{n}{dk}
floor)^2
end{aligned}
]
用 (T=dk) 进行代换，得到：
[egin{aligned}
sum_{d=1}^nd^3sum_{k=1}^{lfloor n/d
floor}mu(k)k^2S(lfloor frac{n}{dk}
floor)^2
&=sum_{d=1}^nd^3sum_{k=1}^{lfloor n/d
floor}mu(k)k^2S(lfloor frac{n}{T}
floor)^2\\
&=sum_{T=1}^nS(lfloor frac nT
floor)^2sum_{d|T}d^3muBig(frac{T}{d}Big)Big(frac{T}{d}Big)^2\\
&=sum_{T=1}^nS(lfloor frac nT
floor)^2T^2sum_{d|T}muBig(frac{T}{d}Big)cdot d\\
end{aligned}
]
根据莫反中的前置知识，有
[varphi(x)=(mu*ID)(x)=sumlimits_{d|x}dcdotmu(frac xd)\\
]
所以：
[=sum_{T=1}^nS(lfloor frac{n}{T}
floor )^2T^2varphi(T)\\
]
然后就可以愉快地杜教筛了。把模板掏下来以便观察：
[g(1)S(n)=sum_{i=1}^nh(i)-sum_{d=2}^ng(d)cdot S(lfloor frac{n}{d}
floor )（杜教筛）\\
(f*g)(x)=sum_{d|x}f(d)g(frac{x}{d})（Dirichlet卷积）\\
varphi*1=ID（和varphi 有关的前置式子）
]
令 (sum(n)=sum_{i=1}^ni^2varphi(i)) ，显然我们需要把 (sum(n)) 中这个 (i^2) 消掉。看着卷积式子，考虑令 (g(n)=n^2) .
于是有：
[(f*g)(n)=sum_{d|n}d^2varphi(d)Big(frac{n}{d}Big)^2=n^2sum_{d|n}varphi(d)
]
然后再把最后一个式子给摁上去：
[(f*g)(n)=n^3.
]
那么就有：
[sum(n)=sum_{i=1}^ni^3-sum_{d=2}^nd^2sum(lfloor frac{n}{d}
floor)\\
]
对于原式：
[sum_{T=1}^nS(lfloor frac{n}{T}
floor )^2T^2varphi(T)\\
]
就可以杜教筛解决了。
//Author: RingweEH
const int N=5e6+10;
int pri[N],tot=0,phi[N],sum[N],Mod;
ll n,inv6,inv2;
bool is[N];
int power( int a,int b )
{
int res=1;
for ( ; b; b>>=1,a=1ll*a*a%Mod )
if ( b&1 ) res=1ll*res*a%Mod;
return res;
}
void Sieve()
{
phi[1]=1; is[1]=1;
for ( int i=2; i<=(N-10); i++ )
{
if ( !is[i] ) pri[++tot]=i,phi[i]=i-1;
for ( int j=1; j<=tot && (i*pri[j]<=(N-10)); j++ )
{
int pos=i*pri[j]; is[pos]=1;
if ( i%pri[j]==0 ) { phi[pos]=1ll*phi[i]*pri[j]%Mod; break; }
else phi[pos]=1ll*phi[i]*phi[pri[j]]%Mod;
}
}
for ( int i=1; i<=(N-10); i++ )
sum[i]=(1ll*i*i%Mod*phi[i]%Mod+sum[i-1])%Mod;
}
int pow1( ll x )
{
x%=Mod; int res=x*(x+1)%Mod*inv2%Mod;
return res;
}
int pow2( ll x )
{
x%=Mod; int res=1ll*x*(x+1)%Mod*(2*x+1)%Mod*inv6%Mod;
return res;
}
int pow3( ll x )
{
int res=pow1(x); res=1ll*res*res%Mod;
return res;
}
unordered_map<ll,int> mp;
int Sieve_Mu( ll x )
{
if ( x<=(N-10) ) return sum[x];
if ( mp[x] ) return mp[x];
int res=pow3(x);
for ( ll l=2,r; l<=x; l=r+1 )
{
r=x/(x/l); res=(res-1ll*(pow2(r)-pow2(l-1))*Sieve_Mu(x/l)%Mod)%Mod;
}
mp[x]=(res+Mod)%Mod; return mp[x];
}
void init()
{
Sieve(); inv2=power( 2ll,Mod-2 ); inv6=power( 6ll,Mod-2 );
}
int main()
{
Mod=read(),n=read(); init();
ll ans=0;
for ( ll l=1,r; l<=n; l=r+1 )
{
r=n/(n/l); ans=(ans+1ll*(Sieve_Mu(r)-Sieve_Mu(l-1))*pow3(n/l)%Mod)%Mod;
}
printf( "%lld
",(ans+Mod)%Mod );
return 0;
}
例题：循环之美
求对于十进制下的 (1leq xleq n,1leq yleq m) ，有多少不相等的小数 (dfrac{x}{y}) 在 (k) 进制下是纯循环小数。
(1leq n,mleq 1e9,kleq 2000) .
Solution
以下过程默认 (gcd(x,y)=1) 。（这个很显然吧）
对于一个 (k) 进制下能成为纯循环小数的数，一定存在一个 (t) ，使得 (k^tcdot dfrac{x}{y}) 的小数部分等于 (dfrac xy) ，即
[k^tcdot dfrac{x}{y}-Biglfloor k^tcdot dfrac{x}{y}Big
floor =dfrac xy =>k^tx-ycdot Biglfloor k^tcdot dfrac xy Big
floor =x\\
]
于是显然有
[k^txequiv x(mod y)=>k^tequiv 1(mod y)\\
k^tequiv 1^t(mod y),k>0,1>0=>kequiv 1(mod y)\\
]
所以有：
[gcd(k,y)=1
]
（这里一开始推到 (y|(k-1)) 去了……我有问题）
这样就可以写出式子了：
[sum_{i=1}^nsum_{j=1}^m[gcd(i,j)=1][gcd(j,k)=1]
]
然后就照常推理：
[egin{aligned}
sum_{i=1}^nsum_{j=1}^m[gcd(i,j)=1][gcd(j,k)=1] &=sum_{j=1}^m[gcd(j,k)=1]sum_{i=1}^n[gcd(i,j)=1]\\
&=sum_{j=1}^m[gcd(j,k)=1]sum_{i=1}^nsum_{d|gcd(i,j)}mu(d)\\
&=sum_{d=1}^nmu(d)sum_{j=1}^{lfloor m/d
floor}[gcd(jd,k)=1]sum_{i=1}^{lfloor n/d
floor}\\
&=sum_{d=1}^nmu(d)Biglfloor frac{n}{d}Big
floor sum_{j=1}^{lfloor m/d
floor}[gcd(jd,k)=1]\\
end{aligned}
]
显然，这个式子应该尽量表示成两个不相关的式子乘积，注意到 (k) 是个常量，于是有
[sum_{d=1}^nmu(d)Biglfloor frac{n}{d}Big
floor [gcd(d,k)=1]sum_{j=1}^{lfloor m/d
floor}[gcd(j,k)=1]
]
令
[F(n)=sum_{d=1}^nmu(d)[gcd(d,k)=1],G(n)=sum_{j=1}^{n}[gcd(j,k)=1]
]
（本来这里的 (F) 是带了一个 $Biglfloor dfrac{n}{d}Big floor $ 的项……但是后来发现这是一个巨大的错误……总之很累赘，所以就没有带，反正最后是要整除分块的，这一项没什么用处）
然后就可以分开考虑了。
[egin{aligned}
F(n)&=sum_{i=1}^nmu(i)[gcd(i,k)=1]\\
&=sum_{i=1}^nmu(i) sum_{d|gcd(i,k)}mu(d)\\
&=sum_{i=1}^nmu(i) sum_{d|k}mu(d)[d|i]\\
&=sum_{d|k}mu(d)sum_{i=1}^{lfloor n/d
floor}mu(id)\\
end{aligned}
]
暴毙了，不会搞了 /kk
……哦，对哦，显然还可以再套一个 ([gcd(i,d)=1]) .不然这个 (mu(id)) 就没了。
这不是更复杂了吗 并不，因为 (mu) 可是积性函数呢 /xyx
[egin{aligned}
F(n)&=sum_{d|k}mu(d)sum_{i=1}^{lfloor n/d
floor}mu(id)\\
&=sum_{d|k}mu(d)sum_{i=1}^{lfloor n/d
floor}mu(i)mu(d)[gcd(i,d)=1]\\
&=sum_{d|k}mu(d)^2sum_{i=1}^{lfloor n/d
floor} mu(i)[gcd(i,d)=1]\\
end{aligned}
]
掉线了掉线了，又不会推了 /kel . 虽然这个 (gcd) 还可以拆但是 (d) 的范围也降不下来了啊。
正在重连……连接失败 . 呜呼！
。去看具体数学了。下午再来推这个。
我回来了！我发现上午人傻掉了。
让我们把最初的样子和最后的结果放在一起，带上一个 (k) 的参数：
[F(n,k)=sum_{d=1}^nmu(d)[gcd(d,k)=1]\\
F(n,k)=sum_{d|k}mu(d)^2sum_{i=1}^{lfloor n/d
floor} mu(i)[gcd(i,d)=1]\\
]
发现下面式子的后半部分就是 (F(lfloor frac{n}{d} floor,d)) ！也就是说可以直接递归计算，然后加个记忆化， (F(n)) 就做完了。
边界为 (F(n,1)) ，杜教筛就好了。筛它就完了!
什么？前面一半漏了？(kleq 2000) 啊，随便跑跑就行了。
下面来看 (G(n)) .
[G(n)=sum_{j=1}^{n}[gcd(j,k)=1]=Biglfloor frac{n}{k}Big
floor varphi(k)
]
没了！这是好的！
那么答案就是：
[Ans=sum_{d=1}^nmu(d)[gcd(d,k)=1]G(lfloor frac md
floor)cdot Biglfloor frac{n}{d}Big
floor
]
整除分块就没了 awa.
注意分子和分母不一样，所以不能 swap(n,m) .
//Author: RingweEH
const int N=5e6+10,K=2010;
int pri[N],tot=0,mu[N],pre_mu[N];
int sum_G[N],n,m,k;
bool is[N];
int gcd( int a,int b )
{
return (b==0) ? a : gcd(b,a%b);
}
void Sieve()
{
is[1]=1; mu[1]=1;
for ( int i=2; i<=(N-10); i++ )
{
if ( !is[i] ) pri[++tot]=i,mu[i]=-1;
for ( int j=1; j<=tot && (i*pri[j]<=(N-10)); j++ )
{
int pos=i*pri[j]; is[pos]=1;
if ( i%pri[j]==0 ) { mu[pos]=0; break; }
mu[pos]=-mu[i];
}
}
for ( int i=1; i<=(N-10); i++ )
pre_mu[i]=pre_mu[i-1]+mu[i];
for ( int i=1; i<=k; i++ )
sum_G[i]=sum_G[i-1]+(gcd(i,k)==1);
}
unordered_map<int,int> sum_Mu,f[K];
int Sieve_Mu( int n )
{
if ( n<=(N-10) ) return pre_mu[n];
if ( sum_Mu[n] ) return sum_Mu[n];
int res=1;
for ( int l=2,r; l<=n; l=r+1 )
{
r=n/(n/l);
res-=(r-l+1)*Sieve_Mu(n/l);
}
return sum_Mu[n]=res;
}
int get_G( int n )
{
return sum_G[k]*(n/k)+sum_G[n%k];
}
int get_F( int n,int k )
{
if ( !n ) return 0;
if ( k==1 ) return f[k][n]=Sieve_Mu(n);
if ( f[k][n] ) return f[k][n];
int res=0;
for ( int i=1; i*i<=k; i++ )
if ( k%i==0 )
{
int x=k/i;
if ( mu[i] ) res+=mu[i]*mu[i]*get_F(n/i,i);
if ( mu[x] && (x*x!=k) ) res+=mu[x]*mu[x]*get_F(n/x,x);
}
return f[k][n]=res;
}
int main()
{
n=read(); m=read(); k=read();
Sieve(); int mn=min( n,m ); ll ans=0;
for ( int l=1,r; l<=mn; l=r+1 )
{
r=min( n/(n/l),m/(m/l) );
ans+=1ll*(n/l)*get_G(m/l)*( get_F(r,k)-get_F(l-1,k) );
}
printf( "%lld
",ans );
return 0;
}