1. 程式人生 > >BZOJ4652: [Noi2016]迴圈之美(莫比烏斯反演,杜教篩)

BZOJ4652: [Noi2016]迴圈之美(莫比烏斯反演,杜教篩)

Description

牛牛是一個熱愛演算法設計的高中生。在他設計的演算法中,常常會使用帶小數的數進行計算。牛牛認為,如果在 k  進位制下,一個數的小數部分是純迴圈的,那麼它就是美的。現在,牛牛想知道:對於已知的十進位制數 n 和 m,在  kk 進位制下,有多少個數值上互不相等的純迴圈小數,可以用分數 xy 表示,其中 1≤x≤n,1≤y≤m,且 x,y是整數 。一個數是純迴圈的,當且僅當其可以寫成以下形式:a.c1˙c2c3…cp-1cp˙其中,a 是一個整數,p≥1;對於 1 ≤i≤p,ci是 kk 進位制下的一位數字。例如,在十進位制下,0.45454545……=0.4˙5˙是純迴圈的,它可以用 5/11
、10/22 等分數表示;在十進位制下,0.1666666……=0.16˙則不是純迴圈的,它可以用 1/6 等分數表示。需要特 別注意的是,我們認為一個整數是純迴圈的,因為它的小數部分可以表示成 0 的迴圈或是 k?1 的迴圈;而一個小 數部分非 0 的有限小數不是純迴圈的。

Input

只有一行,包含三個十進位制數N,M,K意義如題所述,保證 1≤n≤10^9,1≤m≤10^9,2≤k≤2000

Output

一行一個整數,表示滿足條件的美的數的個數。

Sample Input

2 6 10

Sample Output

4
explanation
滿足條件的數分別是:
1/1=1.0000……
1/3=0.3333……
2/1=2.0000……
2/3=0.6666……
1/1 和 2/2 雖然都是純迴圈小數,但因為它們相等,因此只計數一次;同樣,1/3 和 2/6 也只計數一次。

解題思路:

一個喜聞樂見的性質,只要x/y中y與k互質就好了。

所以這道題就是:

$\sum_{i=1}^{N}\sum_{j=1}^{M}\epsilon(gcd(i,j))\epsilon (gcd(j, k))$

$\sum_{j=1}^{M}\epsilon(gcd(j,k))\sum_{i=1}^{N}\epsilon(gcd(i,j))$

$\sum_{j=1}^{M}\epsilon(gcd(j,k))\sum_{i=1}^{N}\sum_{d|gcd(i,j)}\mu(d)$

$\sum_{j=1}^{M}\epsilon(gcd(j,k))\sum_{d=1}^{min(N,M)}\mu(d)\sum_{d|i}^{N}1$

$\sum_{j=1}^{M}\epsilon(gcd(j,k))\sum_{d=1}^{min(N,M)}\mu(d)\left \lfloor \frac{N}{d} \right \rfloor$

${\sum_{d=1}^{min(N,M)}\epsilon(gcd(d,k))\mu(d)\left \lfloor \frac{N}{d} \right \rfloor} \sum_{i=1}^{\left \lfloor \frac{M}{d} \right \rfloor}\epsilon(gcd(i,k))$

莫比烏斯到這裡結束,現在你可以獲得84分,接下來是真正的燒腦環節。

我講的不好,可以看這位巨佬

總之,將後面那個預處理出來。

再二元遞迴求解整體。

程式碼:

  1 #include<map>
  2 #include<cstdio>
  3 #include<algorithm>
  4 typedef long long lnt;
  5 const int N=1000100;
  6 struct pos{lnt x,k;pos(lnt a,lnt b){x=a,k=b;}};
  7 bool operator < (pos a,pos b){if(a.x!=b.x)return a.x<b.x;return a.k<b.k;}
  8 struct Dark_map{
  9     std::map<pos,lnt>A;
 10     void insert(lnt x,lnt k,lnt v){A[pos(x,k)]=v;return ;}
 11     bool hav(lnt x,lnt k){return A.find(pos(x,k))!=A.end();}
 12     lnt val(lnt x,lnt k){return A[pos(x,k)];}
 13 }S;
 14 struct New_map{
 15     std::map<lnt,lnt>A;
 16     lnt a[N];
 17     void insert(lnt p,lnt x){if(p<N)a[p]=x;else A[p]=x;return ;}
 18     bool hav(lnt x){if(x<N)return true;return A.find(x)!=A.end();}
 19     lnt val(lnt x){if(x<N)return a[x];return A[x];}
 20 }Miu;
 21 int prime[N];
 22 int miu[N];
 23 bool vis[N];
 24 int cnt;
 25 int n,m,k;
 26 int twd[N];
 27 int lst[N];
 28 lnt f[2001];
 29 int hd[2001];
 30 lnt gcd(lnt a,lnt b){if(!b)return a;return gcd(b,a%b);}
 31 void adde(int f,int t){cnt++;twd[cnt]=t;lst[cnt]=hd[f];hd[f]=cnt;return ;}
 32 void gtp(void)
 33 {
 34     for(int i=1;i<=k;i++)f[i]=f[i-1]+(gcd(i,k)==1);
 35     for(int i=1;i<=k;i++)for(int j=i;j<=k;j+=i)adde(j,i);
 36     miu[1]=1,cnt=0;
 37     for(int i=2;i<N;i++)
 38     {
 39         if(!vis[i])
 40         {
 41             prime[++cnt]=i;
 42             miu[i]=-1;
 43         }
 44         for(int j=1;j<=cnt&&i*prime[j]<N;j++)
 45         {
 46             vis[i*prime[j]]=true;
 47             if(i%prime[j]==0)
 48             {
 49                 miu[i*prime[j]]=0;
 50                 break;
 51             }
 52             miu[i*prime[j]]=-miu[i];
 53         }
 54     }
 55     for(int i=1;i<N;i++)
 56         Miu.insert(i,Miu.val(i-1)+1ll*miu[i]);
 57     return ;
 58 }
 59 lnt F(lnt x)
 60 {
 61     return (x/k)*f[k]+f[x%k];
 62 }
 63 lnt MIU(lnt x)
 64 {
 65     if(Miu.hav(x))
 66         return Miu.val(x);
 67     lnt tmp=0;
 68     for(int i=2,j;i<=x;i=j+1)
 69     {
 70         j=x/(x/i);
 71         tmp+=1ll*(j-i+1)*MIU(x/i);
 72     }
 73     tmp=1-tmp;
 74     Miu.insert(x,tmp);
 75     return tmp;
 76 }
 77 lnt SUM(lnt Nn,lnt Kk)
 78 {
 79     if(S.hav(Nn,Kk))
 80         return S.val(Nn,Kk);
 81     lnt tmp=0;
 82     if(Nn<1);
 83     else if(Kk==1)
 84         tmp=MIU(Nn);
 85     else{
 86         for(int I=hd[Kk];I;I=lst[I])
 87         {
 88             int x=twd[I];
 89             lnt TMP=miu[x];
 90             if(!TMP)
 91                 continue;
 92             tmp+=SUM(Nn/x,x);
 93         }
 94     }
 95     S.insert(Nn,Kk,tmp);
 96     return tmp;
 97 }
 98 int main()
 99 {
100     scanf("%d%d%d",&n,&m,&k);
101     gtp();
102     lnt ans=0;
103     for(int i=1,j;i<=n&&i<=m;i=j+1)
104     {
105         j=std::min(n/(n/i),m/(m/i));
106         ans+=(SUM(j,k)-SUM(i-1,k))*(lnt)(n/i)*F(m/i);
107     }
108     printf("%lld\n",ans);
109     return 0;
110 }