NOI 2016 循环之美 (莫比乌斯反演+杜教筛)
题目大意:略 洛谷传送门 鉴于洛谷最近总崩,附上良心LOJ链接
任何形容词也不够赞美这一道神题
$\sum\limits_{i=1}^{N}\sum\limits_{j=1}^{M}[gcd(i,j)==1][gcd(j,K)==1]$
$\sum\limits_{j=1}^{M}[gcd(j,K)==1]\sum\limits_{i=1}^{N}[gcd(i,j)==1]$
我们先处理右边的式子$\sum\limits_{i=1}^{N}[gcd(i,j)==1]$:
$\sum\limits_{i=1}^{N}\sum\limits_{d|gcd(i,j)}\mu(d)\sum\limits_{j=1}^{M}[gcd(j,K)==1][d|j]$
$\sum\limits_{d=1}^{N}\mu(d)\left \lfloor \frac{N}{d} \right \rfloor \sum\limits_{j=1}^{M}[gcd(j,K)==1][d|j]$
$\sum\limits_{d=1}^{N}\mu(d)\left \lfloor \frac{N}{d} \right \rfloor \sum\limits_{j=1}^{\left \lfloor \frac{M}{d} \right \rfloor}[gcd(jd,K)==1]$
加下来就是比较关键的一个展开式子,[gcd(jd,K)==1]是[gcd(d,K)==1]&[gcd(j,K)==1]的充分必要条件:
$\sum\limits_{d=1}^{N}[gcd(d,K)==1]\mu(d)\left \lfloor \frac{N}{d} \right \rfloor \sum\limits_{j=1}^{\left \lfloor \frac{M}{d} \right \rfloor}[gcd(j,K)==1]$
84分算法:
令$f(n)=\sum\limits_{i=1}^{n}[gcd(i,K)==1]$
这是啥?
欧拉函数$\varphi$啊!别和我一样反演傻了连欧拉函数的定义都忘了
$f(n)=\left \lfloor \frac{n}{K} \right \rfloor \varphi(K) + \sum\limits_{i=1}^{n\;mod\;K}[gcd(i,K)==1]$
预处理出$\varphi(K)$和$\sum\limits_{i=1}^{n}[gcd(i,K)==1]$,那么$f(n)$就能在$O(1)$求得
可$\sum\limits_{i=1}^{n}[gcd(i,K)==1]\mu(d)$就不能用$\varphi$了,但经过计算,$g$数组能在大约$O(n*K的质因子个数)$的时间内处理出来,极限情况也只不超过$2.6 \cdot 10^{7}$
那么这个问题就在$O(n)$的时间内被解决了
蒟蒻的代码写得非常迷就不放了
100分算法:
这是一道神仙题
先放上原式:$\sum\limits_{d=1}^{N}[gcd(d,K)==1]\mu(d)\left \lfloor \frac{N}{d} \right \rfloor \sum\limits_{j=1}^{\left \lfloor \frac{M}{d} \right \rfloor}[gcd(j,K)==1]$
$\sum\limits_{j=1}^{\left \lfloor \frac{M}{d} \right \rfloor}[gcd(j,K)==1]$可以用$84$分算法里的方法预处理,每次$O(1)$得到
现在令$g(N,K)=\sum\limits_{i=1}^{N}[gcd(i,K)==1]\mu(i)$
$=\sum\limits_{i=1}^{N}\mu(i)*\sum\limits_{d|i}\mu(d)$
$=\sum\limits_{d|K}\mu(d) \sum\limits_{i=1}^{N} [d|i]\mu(i)$
$=\sum\limits_{d|K}\mu(d) \sum\limits_{i=1}^{\left \lfloor \frac{N}{d} \right \rfloor} \mu(id)$
关键的部分来了
如果两个数$gcd(x,y)==1$,说明它们没有公共因子,满足积性函数的性质,那么$\mu(xy)=\mu(x)\mu(y)$
反之它们存在公共因子,那么$\mu(xy)$一定等于$0$
利用这个性质
$=\sum\limits_{d|K}\mu(d)^{2} \sum\limits_{i=1}^{\left \lfloor \frac{N}{d} \right \rfloor} [gcd(i,d)==1]\mu(i)$
诶!右面这个东西$\sum\limits_{i=1}^{\left \lfloor \frac{N}{d} \right \rfloor} [gcd(i,d)==1]\mu(i)$,似乎就是$g(\left \lfloor \frac{N}{d} \right \rfloor,d)$啊!
递归求解即可
当$n==0$是,答案就是$0$
发现$K==1$时不能继续递归了否则会死循环,此时
$G(n,1)=\sum\limits_{i=1}^{n} [gcd(i,1)==1]\mu(i)=\sum\limits_{i=1}^{n} \mu(i)$
$n$可能很大,上杜教筛即可
1 #include <map> 2 #include <cmath> 3 #include <vector> 4 #include <cstdio> 5 #include <cstring> 6 #include <algorithm> 7 #define N1 20000010 8 #define M1 2000010 9 #define K1 2010 10 #define ll long long 11 #define dd double 12 #define cll const long long 13 #define it map<int,int>::iterator 14 using namespace std; 15 16 int N,M,K,maxn; 17 int mu[N1],smu[N1],pr[M1],cnt,phik; bool use[N1]; 18 int f[K1],ps[K1],son[K1],is[K1],pn,sn; 19 vector<int>ss[K1]; 20 void Pre() 21 { 22 int i,j,x; mu[1]=smu[1]=1; 23 for(i=2;i<=maxn;i++) 24 { 25 if(!use[i]){ pr[++cnt]=i; mu[i]=-1;} 26 for(j=1;j<=cnt&&i*pr[j]<=maxn;j++) 27 { 28 use[i*pr[j]]=1; 29 if(i%pr[j]){ mu[i*pr[j]]=-mu[i];} 30 else{ break; } 31 } 32 smu[i]=smu[i-1]+mu[i]; 33 } 34 for(son[++sn]=1,x=K,phik=K,i=2;i<=K;i++) 35 { 36 if(x%i==0){ ps[++pn]=i; phik=phik/i*(i-1); while(x%i==0) x/=i; } 37 if(K%i==0){ son[++sn]=i; } 38 } 39 for(j=1;j<=pn;j++) 40 for(i=ps[j];i<=K;i+=ps[j]) is[i]=1; 41 for(i=1;i<=K;i++) f[i]=f[i-1]+(is[i]?0:1); 42 for(i=1;i<=sn;i++) 43 { 44 x=son[i]; 45 for(j=1;j<=x;j++) 46 if(x%j==0&&mu[j]) ss[x].push_back(j); 47 } 48 } 49 map<int,int>mp; 50 ll Smu(int n) 51 { 52 if(n<=maxn) return smu[n]; 53 it k=mp.find(n); 54 if(k!=mp.end()) return k->second; 55 int i,la;ll ans=1; 56 for(i=2;i<=n;i=la+1) 57 { 58 la=n/(n/i); 59 ans-=1ll*Smu(n/i)*(la-i+1); 60 } 61 mp[n]=ans; 62 return ans; 63 } 64 ll F(int x){return 1ll*(x/K)*phik+f[x%K];} 65 ll G(int n,int k) 66 { 67 if(!n) return 0; 68 if(k==1) 69 { 70 if(n<=maxn) return smu[n]; 71 else return Smu(n); 72 } 73 int i,d;ll ans=0; 74 for(i=0;i<ss[k].size();i++) 75 { 76 d=ss[k][i]; 77 ans+=1ll*G(n/d,d); 78 } 79 return ans; 80 } 81 82 int main() 83 { 84 scanf("%d%d%d",&N,&M,&K); 85 int i,j,la; ll ans=0; maxn=min(max(N,M),20000000); Pre(); 86 for(i=1;i<=min(N,M);i=la+1) 87 { 88 la=min(N/(N/i),M/(M/i)); 89 ans+=1ll*(G(la,K)-G(i-1,K))*(N/i)*F(M/i); 90 } 91 printf("%lld\n",ans); 92 return 0; 93 }