1. 程式人生 > >NOI 2016 循環之美 (莫比烏斯反演+杜教篩)

NOI 2016 循環之美 (莫比烏斯反演+杜教篩)

tar main 說明 pan == amp loj 個數 efi

題目大意:略 洛谷傳送門 鑒於洛谷最近總崩,附上良心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 }

NOI 2016 循環之美 (莫比烏斯反演+杜教篩)