【BZOJ 2671】 2671: Calc （数论，莫比乌斯反演）

2671: Calc
Time Limit: 10 Sec Memory Limit: 128 MB
Submit: 303 Solved: 157
Description

　　给出N，统计满足下面条件的数对(a,b)的个数：
　　1.1<=a<b<=N
　　2.a+b整除a*b

Input

　一行一个数N

Output

　一行一个数表示答案

Sample Input

15

Sample Output

4

HINT

数据规模和约定

Test N Test N

1 <=10 11 <=5*10^7

2 <=50 12 <=10^8

3 <=10^3 13 <=2*10^8

4 <=5*10^3 14 <=3*10^8

5 <=2*10^4 15 <=5*10^8

6 <=2*10^5 16 <=10^9

7 <=2*10^6 17 <=10^9

8 <=10^7 18 <=2^31-1

9 <=2*10^7 19 <=2^31-1

10 <=3*10^7 20 <=2^31-1

Source

【分析】

　　这题的复杂度还挺迷人的。

　　然后$sqrt n$也没发现，以为筛$mu$都要$O(n)$,什么杜教筛的幸好不会。。

　　首先分析$(a+b)|(a*b) → (a/g+b/g)|(a/g*b/g*g) →(a/g+b/g)|g$

　　那就是互质的$a',b'$ 找他们的公倍数$g$就行了。

　　写正常一点就是$$sum_{j=1}^{N}sum_{i=1}^{j-1}dfrac{n}{j*(i+j)} [gcd(i,j)==1]$$

　　到了这里，我就傻眼了，其实嘛。。。j并不会到$n$，只是到$sqrt{n}$

　　$$sum_{j=1}^{sqrt{n}}sum_{i=1}^{j-1}dfrac{n}{j*(i+j)} [gcd(i,j)==1]$$

　　然后我又傻眼了，复杂度迷人的东西啊会把我脑子弄得很乱的。

　　直接枚举j，然后i那里分块，然后就是求[l,r]里面和j互质的数的个数。

　　差分，先求[1,r]里面的，就是$sum_{x=1}^{r}1[gcd(x,j)==1]$

　　即$sum_{d=1}^{r}mu(d)*(r/d)$

　　最后就是$sum_{d=1}^{r}mu(d)*(r/d-(l-1)/d)$

　　枚举约数在前面$sqrt{j}$枚举去了。。

　　真的是暴力出奇迹了。。。

 1 #include<cstdio>
 2 #include<cstdlib>
 3 #include<cstring>
 4 #include<iostream>
 5 #include<algorithm>
 6 #include<cmath>
 7 using namespace std;
 8 #define Maxn 50010
 9 #define LL long long
10 
11 bool vis[Maxn];
12 int pri[Maxn],pl,mu[Maxn];
13 
14 int mymin(int x,int y) {return x<y?x:y;}
15 
16 void init()
17 {
18     memset(vis,0,sizeof(vis));
19     pl=0;mu[1]=1;
20     for(int i=2;i<=Maxn;i++)
21     {
22         if(!vis[i]) pri[++pl]=i,vis[i]=1,mu[i]=-1;
23         for(int j=1;j<=pl;j++)
24         {
25             if(i*pri[j]>Maxn) break;
26             vis[i*pri[j]]=1;
27             if(i%pri[j]==0) {mu[i*pri[j]]=0;break;}
28             mu[i*pri[j]]=-mu[i];
29         }
30     }
31 }
32 
33 int sta[Maxn],sl;
34 
35 void div(int x)
36 {
37     sl=0;
38     int i;
39     for(i=1;i*i<x;i++)
40     {
41         if(x%i==0) sta[++sl]=i,sta[++sl]=x/i;
42     }
43     if(i*i==x) sta[++sl]=i;
44 }
45 
46 int gcd(int a,int b)
47 {
48     if(b==0) return a;
49     return gcd(b,a%b);
50 }
51 
52 int main()
53 {
54     int n;
55     scanf("%d",&n);
56     init();
57     LL ans=0;
58     int sq=(int)(sqrt((double)n));
59     for(int i=1;i<=sq;i++)
60     {
61         div(i);
62         for(int j=1;j<i;)
63         {
64             int x=n/i/(i+j),r;if(x==0) break;
65             r=mymin(i-1,n/x/i-i);
66             for(int k=1;k<=sl;k++)
67             {
68                 ans+=mu[sta[k]]*(r/sta[k]-(j-1)/sta[k])*x;
69             }
70             j=r+1;
71         }
72     }
73     printf("%lld
",ans);
74     return 0;
75 }

View Code

2017-04-06 15:50:26