51nod1222 最小公倍数计数

Posted 2020-09-23 SilverNebula

 

题目来源： Project Euler

基准时间限制：6 秒 空间限制：131072 KB 分值: 640 

定义F(n)表示最小公倍数为n的二元组的数量。

即：如果存在两个数（二元组）X，Y(X <= Y)，它们的最小公倍数为N，则F(n)的计数加1。

例如：F(6) = 5，因为[2,3] [1,6] [2,6] [3,6] [6,6]的最小公倍数等于6。

 

给出一个区间[a,b]，求最小公倍数在这个区间的不同二元组的数量。

例如：a = 4，b = 6。符合条件的二元组包括：

[1,4] [2,4] [4,4] [1,5] [5,5] [2,3] [1,6] [2,6] [3,6] [6,6]，共10组不同的组合。

 

Input

输入数据包括2个数：a, b，中间用空格分隔（1 <= a <= b <= 10^11)。

Output

输出最小公倍数在这个区间的不同二元组的数量。

Input示例

4 6

Output示例

10

 

数学问题 莫比乌斯反演

 

请开始你的反演！
设：

$$ans(n)=\sum_{i=1}^{n} \sum_{j=1}^{n} [\frac{i*j}{gcd(i,j)}<=n]$$
那么 $ans(b)-ans(a-1)$ 就是最终答案

尝试化简上面的式子：
$$\sum_{i=1}^{n} \sum_{j=1}^{n} [\frac{i*j}{gcd(i,j)}<=n]$$
$$\sum_{d=1}^{n} \sum_{i=1}^{\frac{n}{d}} \sum_{j=1}^{\frac{n}{d}} [i*j<=\frac{n}{d}] [gcd(i,j)==1]$$
$$\sum_{d=1}^{n} \sum_{k=1}^{\frac{n}{d}} \mu(k) \sum_{i=1}^{\frac{n}{d}} \sum_{j=1}^{\frac{n}{d}} [i*k*j*k<=\frac{n}{d}] $$
$$\sum_{k=1}^{n} \mu(k) \sum_{d=1}^{\frac{n}{k}} \sum_{i=1}^{\frac{n}{dk}} \sum_{j=1}^{\frac{n}{dk}} [i*j*d<=\frac{n}{k^2}] $$

　　显然d和k值大到一定程度，最后面就是0了，所以我们可以缩小求和上界：

$$\sum_{k=1}^{\sqrt n} \mu(k) \sum_{d=1}^{\frac{n}{k^2}} \sum_{i=1}^{\frac{n}{dk^2}} \sum_{j=1}^{\frac{n}{dk^2}} [i*j*d<=\frac{n}{k^2}] $$

　　这个范围很友好，我们可以枚举$\mu(k)$，求满足条件的i j d三元组数量。
　　需要求的三元组是无序的，为了不重不漏地计数，我们可以分别求出有序（单调上升）的三元组数量，对于其中三个数各不同的、有两个数相同的、三个数都相同的分别计数，然后乘以对应的组合数即可。

 

 1 #include<iostream>
 2 #include<algorithm>
 3 #include<cstring>
 4 #include<cstdio>
 5 #include<cmath>
 6 #define LL long long
 7 using namespace std;
 8 const int mxn=1000010;
 9 int pri[mxn],mu[mxn],cnt=0;
10 bool vis[mxn];
11 void init(){
12     mu[1]=1;
13     for(int i=2;i<mxn;i++){
14         if(!vis[i]){
15             pri[++cnt]=i;
16             mu[i]=-1;
17         }
18         for(int j=1;j<=cnt && pri[j]*i<mxn;j++){
19             vis[pri[j]*i]=1;
20             if(i%pri[j]==0){mu[pri[j]*i]=0;break;}
21             mu[pri[j]*i]=-mu[i];
22         }
23     }
24     return;
25 }
26 LL calc(LL n){
27     if(!n)return 0;
28     LL i,j,k,ed=floor(sqrt(n));
29     LL res=0,tmp=0;
30     for(k=1;k<=ed;k++){
31         if(mu[k]){
32             tmp=0;
33             LL ED=n/(k*k);
34             for(i=1;i*i*i<=ED;i++){
35                 for(j=i+1;j*j*i<=ED;j++)
36                     tmp+=(ED/(i*j)-j)*6+3;
37                 tmp+=(ED/(i*i)-i)*3;
38                 tmp++;
39             }
40             res+=mu[k]*tmp;
41         }
42     }
43     return (res+n)/2;
44 }
45 LL a,b;
46 int main(){
47     init();
48     scanf("%lld%lld",&a,&b);
49     printf("%lld\n",calc(b)-calc(a-1));
50     return 0;
51 }