题解 P1829 【[国家集训队]Crash的数字表格 / JZPTAB】
我的第一篇莫比乌斯反演题解 兴奋兴奋兴奋
贡献一个本人自己想的思路,你从未看到过的船新思路
【分析】
显然,题目要求求的是 \(\displaystyle Ans=\sum_{i=1}^n\sum_{j=1}^m lcm(i,j)\)
根据数论知识,很显然 \(lcm(i,j)={ij\over gcd(i,j)}\)
\(\therefore \displaystyle Ans=\sum_{i=1}^n\sum_{j=1}^m lcm(i,j)\)
\(\displaystyle \qquad\quad=\sum_{i=1}^n\sum_{j=1}^m{ij\over gcd(i,j)}\)
\(\displaystyle \qquad\quad =\sum_{k=1}^{min(n,m)}\sum_{i=1}^n\sum_{j=1}^m{ij\over k}[gcd(i,j)= k]\)
考虑 \(i,j\) 的贡献:
\(\displaystyle \quad Ans=\sum_{k=1}^{min(n,m)}\sum_{i=1}^{\lfloor{n\over k}\rfloor}\sum_{j=1}^{\lfloor{m\over k}\rfloor}{ik\times jk\over k}[gcd(i,j)= 1]\)
由 \(\displaystyle \sum_{d\mid n}\boldsymbol \mu(d)=\sum_{d\mid n}\boldsymbol \mu(d)\boldsymbol I({n\over d})=(\boldsymbol \mu*\boldsymbol I)(n)=\boldsymbol \epsilon(n)=[n=1]\) 得
\(\displaystyle \quad Ans=\sum_{k=1}^{min(n,m)}\sum_{i=1}^{\lfloor{n\over k}\rfloor}\sum_{j=1}^{\lfloor{m\over k}\rfloor}ijk\sum_{d\mid gcd(i,j)}\boldsymbol\mu(d)\)
\(\displaystyle \qquad\quad =\sum_{k=1}^{min(n,m)}\sum_{i=1}^{\lfloor{n\over k}\rfloor}\sum_{j=1}^{\lfloor{m\over k}\rfloor}ijk\sum_{d\mid i\bigwedge d\mid j}\boldsymbol\mu(d)\)
调换顺序,枚举 \(d\)
\(\displaystyle \quad Ans=\sum_{k=1}^{min(n,m)}\sum_{d=1}^{min(\lfloor{n\over k}\rfloor,\lfloor{m\over k}\rfloor)}\boldsymbol \mu(d)\sum_{i=1}^{\lfloor{n\over k}\rfloor}\sum_{j=1}^{\lfloor{m\over k}\rfloor}ijk[d\mid i\bigwedge d\mid j]\)
\(\displaystyle \qquad\quad =\sum_{k=1}^{min(n,m)}\sum_{d=1}^{min(\lfloor{n\over k}\rfloor,\lfloor{m\over k}\rfloor)}\boldsymbol \mu(d)\sum_{i=1}^{\lfloor{n\over k}\rfloor}\sum_{j=1}^{\lfloor{m\over k}\rfloor}ijk[d\mid i][d\mid j]\)
各回各家 ,各找各妈
\(\displaystyle \quad Ans=\sum_{k=1}^{min(n,m)}k\sum_{d=1}^{min(\lfloor{n\over k}\rfloor,\lfloor{m\over k}\rfloor)}\boldsymbol \mu(d)\sum_{i=1}^{\lfloor{n\over k}\rfloor}i[d\mid i]\sum_{j=1}^{\lfloor{m\over k}\rfloor}j[d\mid j]\)
考虑 \(i,j\) 的贡献:
\(\displaystyle \quad Ans=\sum_{k=1}^{min(n,m)}k\sum_{d=1}^{min(\lfloor{n\over k}\rfloor,\lfloor{m\over k}\rfloor)}\boldsymbol \mu(d)\sum_{i=1}^{\lfloor{n\over kd}\rfloor}i\times d[1\mid i]\sum_{j=1}^{\lfloor{m\over kd}\rfloor}j\times d[1\mid j]\)
\(\displaystyle \qquad\quad=\sum_{k=1}^{min(n,m)}k\sum_{d=1}^{min(\lfloor{n\over k}\rfloor,\lfloor{m\over k}\rfloor)}d^2\boldsymbol \mu(d)\sum_{i=1}^{\lfloor{n\over kd}\rfloor}i\sum_{j=1}^{\lfloor{m\over kd}\rfloor}j\)
大部分 dalao 就到此为止,接下来分段用整除分块求解
本人则根据套路又化简了一下:
设 \(\displaystyle T=kd,Sum(n)=\sum_{i=1}^ni={i(i+1)\over 2}\)
\(\displaystyle \quad Ans=\sum_{k=1}^{min(n,m)}k\sum_{T=1}^{min(n,m)}[k\mid T]({T\over k})^2\boldsymbol \mu({T\over k})Sum(\lfloor{n\over T}\rfloor)Sum(\lfloor{m\over T}\rfloor)\)
调换顺序,先枚举 \(T\)
\(\displaystyle \quad Ans=\sum_{T=1}^{min(n,m)}T^2Sum(\lfloor{n\over T}\rfloor)Sum(\lfloor{m\over T}\rfloor)\sum_{k\mid T}{1\over k}\boldsymbol \mu({T\over k})\)
哦吼,凉凉。这是什么鬼东西
一看题目,答案需要对一个大于 \(min(n,m)\) 的质数取模
所以,\({1\over k}\) 显然可以表现为逆元形式
设 \(\boldsymbol {Inv}(n)=n^{-1}\)
\(\displaystyle \quad Ans=\sum_{T=1}^{min(n,m)}T^2Sum(\lfloor{n\over T}\rfloor)Sum(\lfloor{m\over T}\rfloor)\sum_{k\mid T}\boldsymbol{Inv}(k)\boldsymbol \mu({T\over k})\)
大家注意 \(\displaystyle\sum_{k\mid T}\boldsymbol{Inv}(k)\boldsymbol \mu({T\over k})\)
这玩意儿是不是很像迪利克雷卷积呢?我们可以确定 \(\boldsymbol\mu(n)\) 是个积性函数,如果 \(\boldsymbol{Inv}(n)\) 也是积性函数,那它就必然是卷积成一个新的积性函数了
而对 \(\forall gcd(a,b)=1\Leftrightarrow \boldsymbol{Inv}(ab)=(ab)^{-1}\equiv a^{-1}\times b^{-1}=\boldsymbol{Inv}(a)\times \boldsymbol{Inv}(b)(\mod 20101009)\)
所以,我们发现,在取模意义下 \(\boldsymbol {Inv}(n)\) 也是积性函数
所以不妨设 \(\boldsymbol f=\boldsymbol \mu*\boldsymbol {Inv}\)
\(\displaystyle \quad Ans=\sum_{T=1}^{min(n,m)}T^2\boldsymbol f(T)Sum(\lfloor{n\over T}\rfloor)Sum(\lfloor{m\over T}\rfloor)\)
再设 \(\boldsymbol g(n)=n^2\boldsymbol f(n)\)
\(\boldsymbol g(n)=\boldsymbol{id}^2(n)\boldsymbol f(n)\)
所以 \(\boldsymbol g(n)\) 也是积性函数,同时,式子化简完毕:
\(\displaystyle \quad Ans=\sum_{T=1}^{min(n,m)}\boldsymbol g(T)Sum(\lfloor{n\over T}\rfloor)Sum(\lfloor{m\over T}\rfloor)\)
只要我们线性求出 \(\boldsymbol g(T)\) 的前缀和,就可以用整除分块出来了
根据 \(\boldsymbol g(n)\) 是积性函数,我们可以知道它一定是可以线性筛出来的
我们先考虑 \(\boldsymbol f(n)\)
\(\displaystyle \boldsymbol f(1)=\sum_{d\mid 1}\boldsymbol \mu({1\over d})\times d^{-1}=\boldsymbol \mu(1)\times 1^{-1}=1\times 1=1\)
\(\displaystyle \therefore \boldsymbol g(1)=1^2\boldsymbol f(1)=1\times 1=1\)
对 \(\forall p\in prime\)
\(\displaystyle \boldsymbol f(p)=\sum_{d\mid p}\boldsymbol {Inv}(d)\boldsymbol\mu({p\over d})\)
\(\displaystyle \qquad\ =\boldsymbol \mu(p)\times 1^{-1}+\boldsymbol \mu(1)\times p^{-1}\)
\(\displaystyle \qquad\ =(-1)+p^{-1}\)
\(\displaystyle \qquad\ =p^{-1}-1\)
\(\displaystyle \qquad\ =(p^1)^{-1}-(p^0)^{-1}\)
对 \(\forall k>1\)
\(\displaystyle \boldsymbol f(p^k)=\sum_{d\mid p^k}\boldsymbol {Inv}(d)\boldsymbol\mu({p^k\over d})\)
\(\displaystyle\qquad\ =\sum_{t=0}^k\boldsymbol{Inv}(p^{k-t})\boldsymbol \mu(p^t)\)
\(\displaystyle\qquad\ =\sum_{t=0}^1\boldsymbol{Inv}(p^{k-t})\boldsymbol \mu(p^t)+\sum_{t=2}^k\boldsymbol{Inv}(p^{k-t})\boldsymbol \mu(p^t)\)
\(\displaystyle \qquad\ =\boldsymbol{Inv}(p^k)\boldsymbol \mu(1)+\boldsymbol{Inv}(p^{k-1})\boldsymbol \mu(p)+\sum_{t=2}^k\boldsymbol{Inv}(p^{k-t})\times 0\)
\(\displaystyle \qquad\ =(p^k)^{-1}-(p^{k-1})^{-1}\)
\(\therefore \boldsymbol f(p^k)=(p^k)^{-1}-(p^{k-1})^{-1}=(p^k)^{-1}(1-p),k\in N_+\)
\(\therefore \boldsymbol g(p^k)=(p^k)^2\boldsymbol f(p^k)\)
\(\displaystyle \qquad\quad\ \ =p^k\times p^k\times (p^k)^{-1}\times (1-p)\)
\(\displaystyle \qquad\quad\ \ =p^k\times (1-p)\)
\(\displaystyle \qquad\quad\ \ =p^k-p^{k+1}\)
\(\displaystyle \boldsymbol g(p^k)=p\times \boldsymbol g(p^{k-1})\)
\(\therefore\) 对 \(\forall p\in prime\)
\(\boldsymbol g(p\times n)=\begin{cases}p\times \boldsymbol g(n),p\mid n\\\ \\\ \\(1-p)p\times \boldsymbol g(n),p\nmid n\end{cases}\)
即 \(\boldsymbol g(p\times n)=p\times\boldsymbol g(n)\times\begin{cases}1,p\mid n\\\ \\\ \\(1-p),p\nmid n\end{cases}\)
这里对整除分块做一个讲解,已经掌握的 dalao 请跳过:
对于刚刚这个式子,我们希望在已知 \(\boldsymbol g(n)\) 的情况下 \(O(\sqrt n)\) 得出答案
\(\displaystyle \quad Ans=\sum_{T=1}^{min(n,m)}\boldsymbol g(T)Sum(\lfloor{n\over T}\rfloor)Sum(\lfloor{m\over T}\rfloor)\)
这边有一个原理: 对 \(\forall d\in Z_+,\lfloor{n\over d}\rfloor\) 的取值只有不超过 \(2\sqrt n\) 个
证明:对 \(\forall d\leq \sqrt n,{n\over d}\) 只有不超过 \(\sqrt n\) 个取值;对 \(\forall d> \sqrt n,\lfloor{n\over d}\rfloor<\lfloor\sqrt n\rfloor\) ,也只有不超过 \(\sqrt n\) 个取值。加起来,不超过 \(2\sqrt n\) 个
整除分块的原理即从此而来:
我们先考虑单维情况,只含 \(n\) 或只含 \(m\):
设对 \(i\in[l,r]\) 满足:
\(\lfloor{n\over l-1}\rfloor\neq\lfloor{n\over l}\rfloor=\lfloor{n\over i}\rfloor=\lfloor{n\over r}\rfloor\neq \lfloor{n\over r+1}\rfloor\)
则有 \(r=max(i)\Leftrightarrow {n\over r}=min({n\over i})=\lfloor{n\over i}\rfloor\)
\(\therefore r={n\over \lfloor{n\over i}\rfloor}\)
若已知 \(l\) ,则 \(r={n\over \lfloor{n\over l}\rfloor}\)
那 \(l\) 怎么知道呢?
已知第一个为 \(l=1\) ,且每个 \(l\) 对应且仅对应一个 \(r\)
所以当 \(l=r+1\) 时, \(l\) 一定为下个区间的最小值
既然满足如此关系,则:
\(\displaystyle \quad ans=\sum_{T=l}^r\boldsymbol g(T)Sum(\lfloor{n\over T}\rfloor)Sum(\lfloor{m\over T}\rfloor)\)
\(\displaystyle \qquad\quad =\sum_{T=l}^r\boldsymbol g(T)Sum(\lfloor{n\over r}\rfloor)Sum(\lfloor{m\over r}\rfloor)\)
\(\displaystyle \qquad\quad =Sum(\lfloor{n\over r}\rfloor)Sum(\lfloor{m\over r}\rfloor)\sum_{T=l}^r\boldsymbol g(T)\)
设 \(\displaystyle G(n)=\sum_{i=1}^n\boldsymbol g(n)\)
\(\displaystyle ans=Sum(\lfloor{n\over r}\rfloor)Sum(\lfloor{m\over r}\rfloor)[G(r)-G(l-1)]\)
最后只需要 \(Ans=\sum ans\) 进行 \(O(\sqrt n)\) 统计即可
对于二维情况,依葫芦画瓢:
设 \(i\in[l,r],s.t.\)
\(\lfloor{n\over l-1}\rfloor\neq\lfloor{n\over l}\rfloor=\lfloor{n\over i}\rfloor=\lfloor{n\over r}\rfloor\neq \lfloor{n\over r+1}\rfloor\)
\(\lfloor{m\over l-1}\rfloor\neq\lfloor{m\over l}\rfloor=\lfloor{m\over i}\rfloor=\lfloor{m\over r}\rfloor\neq \lfloor{m\over r+1}\rfloor\)
则在 \(l\) 确定的情况下,为了使第一个式子得到满足,会得到 \(r_1={n\over \lfloor{n\over l}\rfloor}\)
为了使第二个式子得到满足会得到 \(r_2={m\over \lfloor{m\over l}\rfloor}\)
故为使得两式同时成立,取 \(r=min(r_1,r_2)\) 即可
剩余的步骤与一维一样 很简单吧
【代码】
那本蒟蒻就放 我码风极丑的 代码了
#include<cstdio>
#include<algorithm>
using namespace std;
#define f(a,b,c,d) for(register int a=b,c=d;a<=c;a++)
#define g(a,b,c,d) for(register int a=b,c=d;a>=c;a--)
//#define LOCAL
typedef int i32;
typedef unsigned int u32;
typedef long long int i64;
typedef unsigned long long int u64;
const i32 MOD=20101009;
const i32 MAXN=1e7+10;
typedef i32 ar[MAXN];
namespace HABIT{//读入输出优化和一堆黑科技,不必管
template<typename T> inline T Max(T a) { return a; }
template<typename T,typename... Args> inline T Max(T a,Args... args){
T b=Max(args...);
return (a>b)?a:b;
}
template<typename T> inline T Min(T a) { return a; }
template<typename T,typename... Args> inline T Min(T a,Args... args){
T b=Min(args...);
return (a<b)?a:b;
}
#ifdef LOCAL
inline char gc() { return getchar(); }
#else
inline char gc() {
static char s[1<<20|1]={0},*p1=s,*p2=s;
return (p1==p2)&&(p2=(p1=s)+fread(s,1,1<<20,stdin),p1==p2)?EOF:*(p1++);
}
#endif
inline i32 read(){
register i32 ans=0;register char c=gc();register bool neg=0;
while(c<48||c>57) neg^=!(c^'-'),c=gc();
while(c>=48&&c<=57) ans=(ans<<3)+(ans<<1)+(c^48),c=gc();
return neg?-ans:ans;
}
char Output_Ans[1<<20|1],*Output_Cur=Output_Ans;
inline bool output() { Output_Cur-=fwrite(Output_Ans,1,Output_Cur-Output_Ans,stdout); }
inline void print(char c) { (Output_Cur-Output_Ans+1>>20)&&output(),*(Output_Cur++)=c; }
inline void print(i32 x){
char buf[20]={0}; (Output_Cur-Output_Ans+sprintf(buf,"%d",x)>>20)&&output();
Output_Cur+=sprintf(Output_Cur,"%d",x);
}
}
using namespace HABIT;
i32 d_N,d_M;
ar ar_d_G,ar_d_Fc,ar_d_Prime;
inline void pre(){
i32 *ptr_d_Prime=ar_d_Prime;
ar_d_G[1]=1;
f(i,2,I,d_N){
if(!ar_d_Fc[i])
ar_d_G[i]=(1ll-i)*(*(ptr_d_Prime++)=(ar_d_Fc[i]=i))%MOD+MOD;
//肯定是负数,取模完直接加上
for(register i32 *p=ar_d_Prime;p<ptr_d_Prime&&*p*i<=d_N;p++){
ar_d_Fc[*p*i]=*p;
ar_d_G[*p*i]=1ll*(*p)*ar_d_G[i]%MOD;
if(*p<ar_d_Fc[i]) ar_d_G[*p*i]=(1ll-*p)*ar_d_G[*p*i]%MOD+MOD;
//肯定是负数,取模完直接加上
else break;
}
ar_d_G[i]+=ar_d_G[i-1];
if(ar_d_G[i]>=MOD) ar_d_G[i]-=MOD;
//前缀和
}
}
int main(){
d_N=read(),d_M=read();
if(d_N>d_M) d_N^=d_M^=d_N^=d_M;
pre();
i32 d_Ans=0;
for(register i32 l=1,r;l<=d_N;l=r+1){
r=Min( d_N/(d_N/l) , d_M/(d_M/l) );
i64 d_Tmp=ar_d_G[r]-ar_d_G[l-1];
i32 d_A=(1ll+d_N/r)*(d_N/r)/2%MOD;
i32 d_B=(1ll+d_M/r)*(d_M/r)/2%MOD;
d_Tmp=d_Tmp*d_A%MOD*d_B%MOD;
d_Ans+=d_Tmp;
if(d_Ans<0) d_Ans+=MOD;
else if(d_Ans>=MOD) d_Ans-=MOD;
}
print(d_Ans);
output();
return 0;
}
最后安利一下 本蒟蒻的博客